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Chair  Philip  L.  Marston 

An  algorithm  to  compute  the  transient  scattering  of  sound  by  elastic  spheres  and 
shells  has  been  developed  and  tested.  Short  tone  bursts  incident  on  elastic  spherical  shells 
under  various  conditions  have  been  investigated.  Previous  experiment  results  for  elastic 
spheres  and  shells  have  been  confirmed.  Some  interesting  physical  aspects  concerning  the 
basics  of  the  scattering  mechanism  have  been  detected.  A  mid-frequency  enhancement  has 
been  examined  for  thin  shells  and  compared  with  the  ray  model  results.  The  ray  model  is 
found  to  closely  approximate  the  enhanced  backscattering  for  tone  bursts  associated  with  a 
subsonic  wave  on  thin  shells.  A  related  filter-like  effect  of  thin  shells  for  the  selected 
frequencies  has  been  detected  and  quantitatively  described  and  the  impulse  response  has 
been  computed.  The  present  computations  can  be  of  practical  use  and  can  be  compared 
with  observations  or  ray  models  of  the  scattering  from  shells  under  the  same  conditions  if 
such  observations  or  models  become  available  later  on.  The  scattering  of  tone  bursts  from 
shells  was  calculated  when  the  carrier  frequency  is  at  or  close  to  the  frequency  of  resonance 
associated  with  a  reverberation  of  longitudinal  waves  across  the  thickness  of  the  shell. 


When  the  shell  is  a  fluid,  the  results  are  easily  interpreted  in  terms  of  an  existing  ray  model. 
When  the  shell  is  elastic,  they  suggest  the  presence  of  a  mechanism  for  prompt  radiation 
noninvolving  complete  circumnavigation  of  the  sphere  by  guided  elastic  waves. 
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CHAPTER  1 
Introduction 

Scattering  from  well-shaped  smooth  objects  such  as  solid  spheres  and  spherical 
shells  has  been  of  interest  and  investigated  for  a  long  dme,  yet  some  basic  physical  aspects 
still  remain  not  well  understood  or  even  unknown.  Recent  attention  has  been  focused  on 
the  geometrical  representation  of  scattering  and  its  extension  to  arbitrary  objects.  For  the 
last  a  couple  of  years.  Dr.  Marston  and  his  research  group  have  done  marvelous  jobs  on  the 
modeling  of  the  scattering  problem  by  ray  synthesis  method  which  has  turned  out  to  work 
excitingly  well.  The  modeling  is,  however,  not  perfect.  There  are  still  some  unclear  physics 
which  presumably  can  be  seen  by  investigating  scattering  of  various  transient  waves. 

There  have  been  numerous  papers  and  researches  on  the  scattering  from  solid 
spheres  and  hollow  spheres,  nevertheless,  relatively  few  transient  experiments  or 
calculations  can  be  found1*6.  One  of  the  difficulties  in  applying  the  geometrical  theory  of 
diffraction  (GTD)  to  elastic  objects  is  that  the  contributions  due  to  surface  guided  elastic 
waves  (SEW),  such  as  Rayleigh  and  Lamb  waves7-8,  can  be  important.  Present  work  is  to 
compute  the  transient  backscattering  from  elastic  spheres  and  spherical  shells  by  using  the 
Fourier  Transform  Method.  This  report  is  written  in  a  consistent  e**401  time  convention 
where  i  =  VT.  Results  verified  the  previous  experiments,  and  most  importantly,  have 
revealed  some  exciting  physics  concerning  the  backscattering  mechanism,  hence  helped 
better  understand  the  scattering  process  and  can  be  helpful  in  further  modeling  of  ray 
method. 
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The  arrangement  of  the  report  is  as  follows.  First,  the  steady  state  scattering  will  be 
briefly  discussed  in  Chapter  2  following  the  introduction.  Chapter  3  presents  the  transient 
scattering  formulation  and  the  method  used  in  the  calculations.  In  Chapter  4,  some  more 
detailed  computational  considerations  can  be  found  and  tests  against  previous 
experimental  work  for  elastic  sphere  will  be  conducted .  Chapter  5  deals  with  the  main  part 
of  the  report  —  backscattering  from  elastic  shells.  Quite  a  few  interesting  results  art 
presented  in  this  chapter.  A  short  conclusion  and  discussion  can  be  found  in  Chapter  6. 
Appendices  will  present  the  programs  used  to  finish  the  job  and  some  illustrations  about  the 
writing  of  the  programs  which  can  be  regarded  as  a  "user’s  guide".  Also  presented  in  the 
appendices  arc  some  additional  data  and  graphs  for  reference  purpose. 
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CHAPTER  2 
Steady-State  Scattering 


2.1  Far  field  scattering 

Figure  1  illustrates  a  generic  problem  of  interest.  The  center  of  an  elastic  sphere  or  a 
smooth  empty  shell  having  a  circular  profile  is  fixed  at  the  origin  of  a  spherical  coordinate 
system. 

Considering  an  infinitely  extended  incident  plane  pressure  wave  of  a  single 
frequency  travelling  in  the  +z  direction 

Pi  (0  *  Po  eite-ta*  *  P0  e  *(z-ct)  ,  (2.1) 

where  Po  is  the  pressure  amplitude ,  k=2 TtlX  is  the  wavenumber,  X  is  the  wavelength  and  c 
is  the  sound  speed  in  the  fluid.  The  exact  PWS  (Partial  Wave  Series)  solution  can  be 
obtained  by  solving  the  differential  equations  subject  to  proper  boundary  conditions.  The 
solution  is  based  on  the  summation  of  infinite  terms  of  spherical  Bessel  and  Hankel 
functions  of  the  first  kind.  This  response  to  a  single  frequency  incident  wave  has  been 
given  a  widely  accepted  name  as  the  reflection  form  function  denoted  here  by  F(ka),  where 
ka  is  the  dimensionless  wavenumber.  Let  F(ka)e’i£0t  be  the  solution  for  constant 
frequency,  then  at  point  P,  which  is  at  a  distance  r»ka2  from  the  center,  harmonic 
scattered  wave  is  produced3 

Pr(t)  =  MM  Po  (2.2) 

'm  i  '  • 


Fig.  1  Geometrical  illustration  of  problem  .  The  incident  pressure  wave 
is  an  infinitely  extended  plane  wave  for  the  steady  state  scattering  and  a 
few  cycles  sine  burst  for  the  transient  scattering.  The  scatterer  can  be 
an  elastic  sphere  or  a  spherical  shell.  The  thickness  of  shell  is  h  »  a  -  b. 
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Eq.(2 2)  has  the  form  of  an  outgoing  spherical  wave  modified  by  the  form  function.  For 
convenience,  we  define  two  dimensionless  variables  as  following: 


x  a  ka  *  (oVc)  a 


T  =  (ct-r)/a  , or  xx  =  k(ct-r). 


(2.3) 


For  the  incident  wave,  r  in  above  expression  should  be  replaced  by  z.  Obviously,  x  is  a 
dimensionless  time  interval,  which  is  analogous  to  retarded  time,  with  a  unit 
comesponding  to  the  time  period  necessary  for  the  wave  to  travel  a  distance  equal  to  the 
radius  of  the  sphere  or  the  out  radius  of  the  shell  in  the  speed  of  sound  in  water. 
Equation(2.2)  can  then  be  rewritten  as 

m = (£)  po  r*  fix) . 


(2.4) 


2.2  Form  Function 


The  complex  scattering  amplitude  or  form  function  in  the  backscattering  direction 
has  the  exact  partial  wave  series  representation 


r*=0 


(2.5) 


where  the  functions  Bn(x)  and  Dn(x)  are  3x3  or  5^5  determinants  obtained  by  satisfying  the 
appropriate  boundary  conditions  for  elastic  sphere9*10, fluid  shell11  or  elastic  shell12. 


The  computation  of  the  form  function  requires  the  truncation  of  the  infinite 
summation  in  Eq.  (2.4).  The  minimum  number  of  terms  (or  the  maximum  number  of  terms 
used  in  the  calculation)  nmax  retained  for  sufficient  convergence  of  F(x)  exceeds  x.  The 
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following  criterion  discussed  by  S.  Kargl  and  Pl.Marston  will  be  used  throughout  this 
work13, 

Umax  *  2  +  [  x  +  4.0x  M3] ,  x<  8 

(2-6> 

nmax  *  3  +  [x  +  ^OSx1^3] ,  x  *  8 
where  the  square  brackets  imply  truncation  to  the  nearest  integer. 
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Chapter  3 


Transient  Analysis 

3.1  Siena!  in  Time  and  Frequency  domain 


For  any  incident  time  signal  g(t),  we  can  write  down  its  frequency  spectrum  as 
Fourier  integral 


and 


g(t)  e+icot  dt 


G((o)e‘i<atd(0 


(3.1) 


(3.2) 


We  can  also  rewrite  these  expressions  in  terms  of  our  dimensionless  variables  x  and  t  as 


GM- 


g(x)  eixT  dt 


(3.3) 


and 


g(T)=27t f  °(X)e'iXXdx 


(3.4) 


I 
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For  each  single  frequency  component  G(x),  we  have  for  its  scattering  response  in 
frequency  domain 

P(x)  *  G(x)  F(x) .  (3.5) 

The  total  far-field  scattering  wave  in  the  time  domain  is  hence  given  by 

P<T)*¥^|_G(x)F(x)rb‘t,h  0.6) 

Since  only  the  backscattering  is  of  interest  in  this  report,  F(x)  in  above  equation  is  the  form 
function  at  specific  angle  0  =*  it. 


then  by  computing  p'(x)  as  a  function  of  x,  the  scattered  wave  as  a  function  of  time 
observed  at  point  P  will  be  obtained.  In  other  word,  a  pulseform  will  be  seen  at  the  point  P 
developing  like  that  computed  according  to  Eq.  3.7  as  value  of  x  increases. 

For  a  real  incident  wave  we  have  g*(x)  =  g(x),  therefore 


G(-x)  =  J  g(x)e*,vtdt 
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J  IgMeidt 


o'M. 


So  we  obtain 


G(-x) »  G*(x) 


Let  the  incident  wave  be  a  Dirac  delta  function  5(t),  then  G(x)  a  1  and 


pMar"  f  F(x)e'l*Tdx 
2tt  J  " 


We  see  that  p’(t)  and  F(x)  are  the  impulse  response  and  the  transfer  function  of  the 
system  (scatterer)  respectively.  They  are  a  Fourier  Transform  pair.  Since  p’(t)  is  a  real 
function,  following  the  above  derivation,  we  immediately  write  down 

F(-x)  -  F*(x)  (3.9) 

Now,  we  are  ready  to  simplify  our  Eq.(3.7)  as 


PW-  2jj  I  G(x)F(x)e-ixTdx 


»2 „f  C3(x)F(x)e*ixXdx+^|  G(x)  F(x)e-ixt  dx 
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--L 

lit 


j  G(-x)n-x)CixTdx+^^ 


G(x)F(x)rixTdx 


G(x)F(x)e-ixxdx  +  -L 


G(x)  F(x)e-ixtdx 


G(x)F(x)e-ixTdx 


(3.10) 


Eq.  (3.10)  is  the  basic  formula  we  will  use  throughout  the  work.  By  using  Eq.(3.10),  we 
don’t  need  to  worry  about  the  physical  meaning  of  negative  frequency.  This  equation  has 
turned  out  to  be  important  in  understanding  the  algorithm  developed  later  by  using  Digital 
Fourier  Transform  method. 
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where  coo  is  the  carrying  frequency  of  the  burst  and  xo  =  koa  =  (Ooa/c  is  its  dimensionless 
counterpart  Its  spectrum  can  be  calculated  as  follows. 
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GW-  JC  g(t)  c‘*t  dt 


where  a  =  x+xo,  P  =  x-xo-  The  integration  in  Eq.(3.10)  was  accomplished  by  utilizing 
the  simple  trapezoidal  method.  To  assure  the  detail  of  the  scattered  wave  be  properly 
displayed,  the  interval  between  successive  t  values  has  to  be  chosen  small  enough.  This 
value  has  been  set  to  At=  so  that  there  are  25  points  over  one  cycle  of  sine  bunt  (see 
more  detail  in  Appendix  A6).  For  the  case  of  a  4-cycle  bunt  with  x0  =■  49.1,  At »  0.02 
and  xmax  was  chosen  to  be  400.  This  value  of  xmax  is  determined  by  the  fact  that  the 
inverse  transformed  vereion  has  been  good  enough  for  a  cut-off  frequency  of  400  and  the 
convergence  beyond  400  is  so  slow  that  in  order  to  get  significant  improvement  over  the 
case  of  xmax  »  400,  xmax  has  to  be  greater  than  1000,  which  implies  whole  lot  a  more 
computational  work.  The  ideal  incident  tone  bunt  and  computed  backscattering  from  an 
object  with  F(x)  *  1,  or  the  real  incident  burst  are  both  shown  in  Fig.  2. 

Replacing  the  integral's  upper  limit  with  a  finite  value  is  the  same  as  letting  the 
signal  pass  through  an  ideal  low-pass  filter.  Therefore  a  filtered  replica  of  the  ideal  incident 
burst  is  expected.  This  effect  is  clearly  shown  by  the  actual  incident  burst  depicted  in  Fig.  2 
where  the  sharp  part  of  the  ideal  burst  has  been  rounded  off  since  the  high  frequency 
information  above  the  xmax  have  been  lost.  But  the  amplitude  of  those  high  frequency 
components  are  small  enough  that  the  difference  is  not  significant  as  far  as  the  transient 
scattering  is  concerned.  Physical  aspects  are  still  there  and  almost  unaffected. 
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42  Test  against  wneriment 

The  test  against  the  experiment  result  obtained  by  William  and  Marston9  for  a 
tungsten  carbide  sphere.in  water  was  then  carried  out.  The  material  parameters  used  for 
tungsten  carbide  were:  density  p*  13.80  g/cm^,  longitudinal  wave  speed  q=  6.860  *10  5 
cm/s,  shear  wave  speed  =  4.185  *  10  5  cm/s.  Those  for  water  were:  density  pw=  lg/cm3, 
longitudinal  wave  speed  =  1.4760  *  10  ^  cm/s.  The  incident  burst  is  again  chosen  as  stated 
in  Eq.  42  with  b=4.  The  form  function  F(x),  calculated  on  the  Dec  station  in  the  range  of 
0<x<400  with  an  interval  Ax  =  0.05  and  plotted  on  HP-7475A  plotter,  is  shewn  as  Fig.3b 
for  completeness.  Both  the  experiment  and  calculated  results  are  also  shown  in  Fig.  3 
where  we  have  seen  that  they  match  each  other  reasonably  well.  The  form  function  is 
referenced  at  the  center  of  the  sphere  so  that  scattered  wave  leads  the  incident  wave  by  t 
=  2.  This  test  simply  verified  that  the  approach  to  the  transient  problem  is  successful. 

4.3  DFT  .Approach 

To  do  the  computation  more  effectively  we  utilize  the  FFT  programs  provided  by 
the  CMS  main  frame.  To  better  understand  the  algorithm  we  are  using  we  have  to  talk 
about  the  DFT  (abbreviation  for  Digital  Fourier  Transform)  first 

DFT  14  of  a  sequence  of  N  samples,  f(nT),  0<n<N-l,  is  defined  as  another 
sequence 

mm 

F(jQ)=  £  f(nT)eiflTni 

n  =0 


(4.4) 
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Fig. 3a.  Experiment  result  of  the  backscattering  from  a  tungsten  carbide  sphere  in 
water.  The  incident  4-cycle  sine  burst  has  a  carrier  frequency  xo  =  49.1.  The 
individual  echoes  are  labeled  as  follows:  A-  specular  reflection,  B-Rayleigh 
surface  wave  echoes,  C-other  surface  waves  and  transmitted  bulk  waves. 
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The  inverse  DFT  is  given  by 


W-iZ  J<jQ)c-iOT"j 

J-0 

where  0£n£N-l,0£j£N-l  and 


NT 


(4.5) 


(4.6) 


Clear  enough,  N  has  to  be  very  large  if  we  intent  to  have  both  frequency  and  time 
intervals  to  be  sufficiently  small  so  that  we  don’t  lost  any  detailed  information  during  the 
computations. 


4.4  Sampling 

Now  we  consider  an  input  sequence  Yj  obtained  by  sampling  the  incident  sine  bum 
at  an  equal  interval  Ax.  Since  DFT  regards  the  input  sequence  as  one  period  of  a  periodic 
signal,  we  let  Yj  be  described  as 

sin  (x0Bm/M)  m=0,l, . M-l 

Yj=  {  (4.7) 

0  m=*  M,  M+1.....N , 

where  M  is  the  number  of  points  over  the  duration  of  the  burst  while  N  is  the  length  of  the 
input  sequence.  A  look  at  the  measurements  of  backscattering  by  elastic  sphere  and  shells 
tells  us  that  the  ratio  of  N/M  must  be  greater  than  16  if  we  want  to  distinguish  at  least  the 
first  few  scattered  surface  waves.  By  letting  the  ratio  big  enough,  we  can  assure  die  period 
of  the  input  sequence  to  be  long  enough  to  cover  the  transient  time  interval  we  are  interested 
in.  Otherwise,  the  second  period  of  the  specular  reflection  would  get  back  to  our 
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measurement  point  P  before  the  SEW  due  to  the  first  burst.  This  is  the  intrinsic 
disadvantage  of  the  method  over  the  previous  one  where  we  only  have  one  period  of  burst 
This  turns  out  not  to  be  a  problem  because  in  order  to  satisfy  other  conditions  to  be 
discussed  below,  N  has  to  be  very  large  compared  with  M,  which  ensures  the  above 
condition  in  any  case. 

Sampling  period  At  has  to  be  chosen  so  that  there  are  sufficient  sampling  points 
over  one  cycle  of  the  burst  And  at  the  same  time  we  need  to  be  sure  that  the  frequency 
interval  Ax  is  small  enough.  But  calculating  of  the  form  function  involves  a  lot  of 
computing  time.  Ax  is  therefore  set  to  be  0.05  which  has  been  proven  sufficient  by  the 
later  computational  results.  Actually  one  can  not  ordinarily  discern  the  results  of  Ax  =  0.05 
and  Ax  =  0.01 .  Suppose  we  want  L  points  per  cycle,  then, 

M  =  LX  b 

At  *  -22 L. 

L*o 

and 

At  Ax  Ax  (4.10) 

for  L  »  20  and  xo=49.1,  N=  19,640.  The  maximum  x<j  value  involved  in  this  work  is 
about  340  and  resulting  huge  number  of  N  can  be  a  computing  problem.  To  avoid  this  we 
will  have  to  decrease  the  value  of  L  and  increase  the  value  of  Ax  sometimes  ( see  more 
detail  in  the  Appendix  A6 ). 


(4.8) 

(4.9) 
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4.5  FFT  Algorithm 

FFT  programs  we  used  are  FFTCF  and  FFTCB  routines  in  the  IMSL  bank  1 5. 
FFTCF  routine  computes  the  discrete  complex  Fourier  Transform  of  a  complex  vector  of 
size  N  while  FFTCB  computes  the  inverse  discrete  transform  of  a  complex  vector  of  size 
N.  The  method  used  is  a  variant  of  the  Cooley-  Tukey  algorithm  which  is  most  effective 
when  N  is  a  product  of  small  prime  factors.  Specifically,  given  an  N-vector  Y,  FFTCF 
returns  another  N-vector  Z 


N-l 

Zj-£  Y„  e2R“*i/N ,  j  =  0,1,...  N-l 

n»0 

and  given  an  N-vector  Z,  FFTCB  returns  another  N-vector  Y 

n-i 

Yn  -  jJ-X  ZjC-^VN.  n  *  0,1,...  N-l 
j-0 

For  a  real  sequence  Yn*  =  Yn,  we  have 
N-l 

ZN.j  =  X  Yn  e^N-jVN 

n  *  0 

Yn  e*2*ini/N 

n  «*  0 
N-l 

-  Z  (Yn  e*N/Nj* 

n  ■  0 


(4.11) 


(4.12) 


or. 
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Ztf.f*  Zj*  (4.13) 

Above  equation  is  analogous  to  Eq.(3.8).  Clearly,  Zj  ( j>[N/2J)  corresponds  to  the  negative 
frequency  components  and  is  the  highest  frequency  component  calculated  by  the  FFT. 
Therefore  we  can  also  write 

n=»0,l,...,N-l  (4.14) 


n  =  0.1....JN-1  (4.15) 

The  problem  of  reference  point  can  be  serious  for  the  FFT  calculation.  If  we  still  let 
the  incident  burst  be  described  by  Eq.(4.7),  we  would  expect  the  backscattering  echoes 
leading  the  incident  burst  by  x  =  2.  In  other  word,  backscattering  would  start  at  x  =  -2.  But 
FFT  routines  do  not  provide  any  data  in  the  negative  time  scale.  In  fact,  since  the  incident 
burst  is  regarded  as  one  period  of  a  periodic  input  sequence,  the  echo  starting  at  x  =  -2 
would  appear  at  the  end  of  the  first  period  of  computed  backscattering.  Basically,  there  are 
two  ways  to  solve  the  problem.  One  way  is  to  multiply  each  frequency  component  of 
incident  burst  by  a  phase  shift  factor  exp(i2ka)  before  doing  the  integration  or  summation. 
This  is  equivalent  to  taking  the  convolution  of  the  incident  burst  with  S(x-2).  But  this  8 
function  could  not  be  perfect  due  to  the  upper  limit  of  integration  or  summation.  This  will 
cause  some  distortion  on  the  real  incident  burst.  The  other  way  is  to  directly  shift  the 
incident  burst  over  the  time  scale  so  that  is  starts  at  x  =  2  instead  of  x  =  0.  This  is  in 


(  [N/2]  \ 

fn=2Rd-l.  £  Zj 
\  j«o 


Eq.(3.10)  becomes 


/  [N/2] 

3'„*2Rd£  2  ZjFj  e-a^i/N 
j-o 
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principle  the  same  as  the  first  method  but  avoids  the  other  unwanted  effect  In  order  to  see 
the  early  stage  of  the  backscattering  clearly,  we  have  actually  shifted  the  incident  burst  so 
that  it  starts  at  t  ■  3.  The  backscattering  echoes  therefore  starts  at  t  *  1  in  any  case.  The 
ideal  and  real  (inverse  transformed  assuming  a  hypothetical  body  with  F(x)  m  1)  incident 
bunts  for  FFT  routines  are  shown  in  Hg.4a.  The  spectrum  of  the  burst  was  computed  by 
P'FTCF  routine  and  the  backscattered  wave  was  computed  by  FFTCB  routine  according  to 
Eq.(4.15)  where  we  have  let  Zj  ■  0  for  j  *  [n/2]+l, ...  N-l.  As  we  see  from  Fig.4b,  The 
backscattered  wave  calculated  by  using  the  OFT  is  exactly  the  same  as  our  previous 
computed  result  by  doing  the  integration  directly  (Fig.3). 


! 


OH 
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Chapter  5 


Backscattering  From  Elastic  Shells 
5.1  Verification  with  Experiments 

The  Sana  FFT  algorithm  was  applied  to  the  case  of  a  16.2%  -  thick  spherical  shell. 
The  material  constants  are  listed  in  Table  2.  The  results  for  ka  »  68.8,  64.7,  36.4  are 
presented  in  Fig.  5.  Incident  burst  was  again  chosen  to  be  four-cycle  sine  wave.  These 
results  reasonably  match  the  previous  experiment  result  obtained  by  S.  Kargl  and  P.L. 
Marston12  which  is  also  shown  in  Fig.5  for  comparison.  Evidently,  there  are  some 
significant  difference.  In  the  computation  case  a  slow  varying  wave  packet  shows  up, 
which  was  not  expected  and  understood  at  the  very  beginning  of  the  report  at  all.  But  we 
can  get  some  clue  by  looking  at  the  form  function  of  such  a  shell  which  is  also  given  in 
Fig.5.  Clearly  there  is  a  big  hump  at  the  very  low  frequency.  This  phenomenon  explains,  at 
least  tentatively,  why  the  low  frequency  components  of  the  incident  wave  were  sort  of 
amplified  and  the  effect  was  so  evident.  The  central  dimensionless  wavenumber  of  the 
hump  is  around  x  *  5.  For  the  a  *  1.905cm  shell,  it  corresponds  a  low  frequency  of  fh  = 
62kHz  which  is  much  lower  than  the  actual  carrier  frequency  of  the  incident  burst.  This 
enhancement  of  the  low  frequency  response  may  have  been  suppressed  during  the 
experiments  since  the  combined  transmitting  and  receiving  apparatus  would  not  have  had 
the  broad  bandwidth  assumed  in  the  calculation.  It  is  especially  plausible  that  the  four-cycle 
source  burst  used  in  the  experiments  did  not  have  the  low  frequency  spectral  components 
implicit  in  Eq.(4.3).  For  thinner  shells,  the  enhancement  is  shifted  to  a  higher  frequency 
and  becomes  what  is  known  as  the  mid-frequency  enhancement  of  the  backscattering. 


Material 


p  (g/cm^)  q  (km/s)  cs  (km/s) 


440  C  stainless  stell  7.84  5.854  3.150 

Water  1.00  1.479 

rssaaaaasagaasassaggaaarasaaaasaasaasassagass; 

Table  2.  Material  parameters  for  SS440C  shell  and  water. 


FORM  FUNCTION  FOR  SPHERICAL  SHELL 


FORM  FUNCTION  FOR  SPHERICAL  SHELL 
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52  More  on  Mid-freouencv  Enhancement 


The  unexpected  phenomenon  and  the  tentative  yet  reasonable  explanation  lead  us  to 
the  further  examination  of  some  other  shell  cases.  Because  the  hump  effect  will  most  likely 
happen  for  thin  shells,  the  backscattering  of  a  4-cycle  sine  burst  incident  on  a  5%  -  thick 
spherical  shell  was  computed.  The  carrier  frequency  Xq  =  66.3  was  chosen,  which  happens 
to  be  a  form  function  minimum.  Shown  in  Fig.6  are  the  form  function,  spectrum  of  the 
incident  burst  and  the  transient  backscattering.  This  time  we  see  a  huge  hump  or  "hill”  in 
the  range  of  ka  values  18-24  in  Fig.6a  and  a  relatively  large  amplitude  of  components  in 
this  range  for  the  incident  burst  in  Fig.6b.  These  combined  factors  clearly  explain  the  low 
frequency  signal  detected  in  Fig.6c.  There  is  a  relatively  large  response  to  the  initial  weak 
low  frequency  component  of  the  incident  bunt 

In  the  system  point  of  view,  the  shell  is  acting  like  a  filter  with  a  big  response  for 
ka  =  18-24.  Any  significant  frequency  component  in  this  range  for  any  input  signal  will  be 
singled  out  Then  what  about  a  bunt  with  a  carrier  frequency  in  that  range?  The  question 
is  well  answered  by  the  huge  echo  following  the  specular  reflection  as  shown  in  Fig.6d  for 
a  4-cycle  sine  bunt  with  xo  *  22. 1.  This  particular  value  of  x<j  was  picked  up  because  in 
Fig.6c  we  see  the  low  frequency  response  looks  like  a  sine  wave  of  1/3  of  original  sine 
signal  carrier  frequency.  The  existence  of  a  mid-frequency  enhancement  or  "hill”  in  the 
form  function  for  the  steady  state  backscattering  of  elastic  shells  has  been  noted  by  various 
authors  including  Taimant,  Uberall,  Miller,  and  Dickey*  and  Sammelmann,  Trivett,  and 
Hackman16*17.  The  effect  is  interesting  enough  that  it  deserves  further  investigation. 
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53  Mid-freauencv  Enhancement:  Rav  model  and  Verification 

First  of  all,  we  notice  the  arrival  time  of  the  echo.  Relative  to  the  specular  response, 
this  particular  echo  always  has  a  time  delay  of  Xd  ■  4-5,  that  is  the  time  period  necessary  for 

the  signal  to  travel  a  distant  equal  to  the  4  or  5  times  of  out  radius  of  the  shell  in  terms  of 
the  speed  of  sound  in  water.  To  facilitate  the  understanding  of  the  major  scattering 
mechanisms,  the  ray  model  analysis  borrowed  from  Marston  .Kargl  and  Williams11*18*20 
is  simply  restated  here.  Fig.7a  shows  the  ray  diagram  for  leaky  Lamb  waves  where  a  trace- 
velocity  matching  condition  gives  for  c/  >  c 

0/  *  arcsin  (c/c/) .  (5.1a) 

Ray  synthesis  of  form  function  has  basically  two  parts,  one  is  due  to  the  specular 
reflection  and  the  other  is  the  contribution  of  the  variety  of  SEW.  It  can  be  written  in  the 
following  form  where  /  denoted  the  class  of  SEW: 

f  =  fsp  +  f/=s0  +  U- ao  +f/=ao' +  •••  *  (5.1b) 

where  ff=ao'  term  is  thought  to  be  especially  important  for  thin  shells.  Some  authors16*17 
use  to  denote  the  wave  ao'.  The  geometric  series  for  repeated  leaky  wave  contributions 
gives  for  the  steady  state  scattering  from  a  spherical  shell18*19 
f  _  -Gjcxp(-2(tt  -  0/)  P<]  exp (iq/) 

1  +exp(-2rtP/)exp(i2jtkac/c/)  >  (5-2) 


where 


T|/=  2ka[(c/c/)  (Jt  -  0/)  -  cos0/]  -  tc/2 
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Fig.7a.  Ray  diagram  for  contributions  to  backscattering  due  to  various  kinds  of 
surface  wave  ( of  type  0  excited  on  an  elastic  spherical  shell.  Ray  path  ABECD  is 
the  usual  case  where  c/  >  c.  Ray  path  A'B'E  CD’  and  AA"B"E  C"D"  correspond 
to  ci  =  c  (known  as  "creeping  wave")  and  q  <  c  (known  as  "trapped  wave”).  In  the 
esae  of  c fi  <  0,  The  energy  follows  the  short  path  of  ABCD  without 
circumnavigating  the  shell. 
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Fig. 7b.  The  ray  diagram  for  the  specular  reflection  from  an  elastic  spherical  shell. 
The  point  S  is  the  vertex  of  the  refracting  surface.  The  ray  AB  infinitesimally  close 
to  the  z-axis  is  incident  on  the  shell  with  impact  parameter  s.  The  ray  ABC  is  the 
ordinary  specular  ray:  while  ABDEF  and  ABDEGHI  are  the  first  two  internal 
specular  reflection  contributions.  Intersection  of  the  projection  (dashed  lines)  of  the 
outgoing  rays  and  the  z-axis  define  locations  of  virtual  point  sources,  Vn,  which 
describe  local  curvature  of  the  wavefront  associated  with  each  ray.  The  specular 
reflection  contribution  to  the  form  function  for  backscattering  is  determined  from  a 
superposition  of  the  wavefields  from  the  virtual  sources. 
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IGil «  8«P/c/c/  .  (5.3) 

Some  comments  on  the  nature  of  the  specular  reflection  are  appropriate.  Fig.7b 
shows  the  ray  diagram  for  the  external  specular  reflection  and  the  other  specular 
contributions  associated  with  the  reverberation  of  longitudinal  waves  within  the  shell.  The 
point  S  is  the  vertex  of  the  refracting  surface.  The  ray  AB  is  infinitesimally  close  to  the  z- 
axis.  The  ray  ABC  is  the  external  specular  ray;  while  ABDEF  and  ABDEGHI  are  the  first 
two  internal  specular  reflection  contributions.  At  the  internal  surface,  the  rays  are  totally 
reflected  and  the  rays  are  partially  reflected  and/or  transmitted  at  the  water-shell  interface. 
Intersection  of  the  projection  (dashed  lines)  of  the  outgoing  rays  and  the  z-axis  define 
locations  of  virtual  point  sources,  Vn,  which  describe  local  curvature  of  the  wavefront 
associated  with  each  ray.  For  thin  shells  numerical  computations  of  steady-state 
scattering19  suggest  that  it  is  satisfactory  to  neglect  the  multiple  reverberations  in  the  region 
of  the  mid-frequency  enhancement  and  to  approximate  the  specular  contribution  is  that  due 
to  ray  ABC 

Fig.8,9,10  calculated  by  Dr.  Sun  for  a  2.5%-thick  SS304  shell  show  how  the 
radiation  damping  parameter  (3/  (in  np/radian),  the  phase  velocity  c /  and  the  group  velocity 
Cg i  of  ao,so  and  a©  waves  vary  with  ka.  For  the  thickness  of  the  shell  we  arc  considering 
here  these  are  the  only  relevant  surface  guided  waves.  The  method  of  computing  these 
surface  wave  properties  is  based  on  applying  the  Watson  transform  methodology  to  the 
exact  elastic  equations  for  shells  as  discussed  in  Refs.  12,13, and  20.  In  Fig.8,  we  see  that 
in  the  region  of  34<x<62,  P/  is  near  zero  for  ao  and  is  large  for  s0  wave.  Eq.(5.2)  and 
(5.3)  tell  us  that  the  contribution  to  the  form  function  due  to  a  surface  wave  will  be  small  if 


short-long-dashed  curve  is  for  /  = 
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P<  value  of  the  surface  wave  is  either  small  (coupling  coefficient  G /  will  be  small)  or  large 
(big  damping  will  kill  it),  hence  only  the  ao'  wave  is  important  in  this  analysis.  Fig.9  tells 
us  that  the  phase  velocity  for  ao’  in  the  ka  range  we  are  dealing  with  is  very  close  to  the 
sound  speed  in  water  c,  which  means  that  we  are  near  the  coincident  condition.  Eq.(S.lb) 
immediately  gives  0/  =  n/2  for  c/  =  c.  In  Fig.7a,  a  traveling  path  of  the  wave  is  drawn  as 
A'B'E  CD*.  The  wave  first  travels  with  speed  of  c  straight  along  the  ray  A'B'  and 
couples  with  the  surface  of  the  shell  in  the  vicinity  of  point  B',  then  goes  with  speed  of  Cg / 
along  a  ray  path  BE  C  that  lies  along  the  periphery  of  the  shell.  Finally  travels  to  D’  with 
speed  c  straight  along  the  path  CD'.  This  is  the  case  of  "creeping  wave".  For  c/<  c,  a 
"trapped  wave"  exist  where  the  wave  is  tunnelling  into  the  shell  in  the  form  of  an 
evanescent  wave  as  indicated  by  the  rays  A"B"E  CD"  in  Fig. 7a.  The  terminology  is  that  of 
Ho  and  Felsen21.  For  both  cases,  the  incident  angle  is  taken  to  be  n/2.  The  applicability  of 
Eq.(5.1b)  to  steady  state  scattering  from  thin  shells  was  recently  demonstrated  by  Sun  and 
Marston19.  Inspection  of  Fig.9  shows  that  c/  <  c  for  the  ao'  wave  which  exhibits  "trapped" 

wave  behavior  at  low  frequencies  becoming  "creeping  wave"  like  as  ka  is  increased  above 
-40. 


Consider  now  the  implication  of  Fig.7a  to  the  scattering  of  tone  bursts.  The  arrival 
time  of  the  burst  follows  from  the  inspection  of  Fig.7a.  Since  9/  =  7t/2  for  creeping  and 
trapped  waves  the  time  delay  of  ao'  wave  t<j  relative  to  the  specular  reflection  in  terms  of 
the  dimensionless  time  interval  can  be  approximated  as 

Td  =  (l+Jtc/cgf+l),  (5.4) 

where  the  actual  time  delay  is  td  =  x^djc.  Fig.  10  states  that  the  group  velocity  for  ao'  is 
close  to  and  a  little  bit  bigger  than  c.  In  the  region  of  34<x<62,  eg i/c  »  1.38-1.03,  which 
corresponds  to  a  t<j  value  of  4.28  -  5.05.  The  property  of  c /  and  eg /  for  ao’  wave  was 
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found  to  be  common  to  other  thin  shells  and  hence  the  above  analysis  should  also  apply  to 
the  general  arrival  time  of  the  ao'  wave  bunt  for  the  5%-thick  SS440c  shell.  We  see  that 
this  simple  approximation  agrees  with  our  previous  observation  in  Fig.6. 

To  further  test  the  validity  of  this  ray  model  and  investigate  die  backscattering  under 
the  mid-frequency  enhancement  condition,  we  have  made  the  computations  for  the 
backscattering  of  bunts  with  various  carrier  frequencies.  Since  the  available  ray  method 
calculation  by  Sun  is  for  the  case  of  a  2.5%  -  thick  SS304  shell,  the  same  parameten  were 
chosen  and  listed  in  table  3.  The  carrier  frequencies  picked  up  were  ka  *34,  38,  42,  46, 
50,  54,  58,  62.  In  order  to  measure  the  correct  amplitude  for  the  transient  waves  the 
incident  bunt  has  to  be  long  enough  so  that  the  near  steady  state  response  can  be  achieved 
before  the  signal  turns  off.  On  the  other  hand,  it  could  not  be  too  long,  otherwise 
additional  echoes  would  come  in  and  make  the  measurement  impossible.  The  establishment 
of  the  steady  state  response  has  turned  out  to  be  slow;  as  a  compromise,  a  20-cycle  bunt 
was  finally  chosen.  The  calculated  transient  backscattering  arc  depicted  in  Fig.  11.  For  each 
case,  the  amplitude  of  this  particular  type  of  wave  relative  to  the  specular  echo  were 
measured  by  looking  at  the  data  presented  in  Table  4  and  averaging  the  closest  four 
maximum  values  highlighted  by  the  rectangular  boxes.  Table  5  provides  the  predicted 
amplitude  value  by  the  ray  mode!  and  that  of  present  calculation  for  each  ka  value  we 
picked  up.  The  amplitude  measured  from  the  transient  backscattering  is  depicted  by  the 
stan  in  Fig.  12  where  the  ray  model  calculation  has  been  drawn  as  a  solid  line  in  the  same 
figure  for  comparison.  The  ray  model  calculation  was  based  on  the  following  equation 

|f/»«wn»o|  *  8rtP^exp(-P/jr)  (5  5) 

where  m  *  0  indicates  that  the  contribution  is  only  due  to  the  first  radiation  after 
circumnavigating  half  the  circle  which  is  suitable  to  compare  with  the  transient 


Material 


p  (g/cm^)  ci  (km/s)  cs  (km/s) 


SS304 

Water 


7.570  5.675  3.141 

1.00  1.479 


Table  3.  Material  parameters  for  stainless  steel  304  shell  and  water. 
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T 

P 

T 

P 

6.596870 

0.036035 

6.606110 

0.111537 

6.615349 

0.159312 

6.624588 

0.209079 

6.633827 

0.235192 

6.643066 

0.232162 

6.652306 

0.223170 

6.661545 

0.178535 

6.670784 

0.126900 

6.680023 

0.062583 

6.689262 

-0.012010 

6.698502 

-0.075965 

6.707741 

-0.142101 

6.716980 

-0.184319 

6.726219 

-0.216449 

6.735458 

-0.222948 

6.744698 

-0.208140 

6.753937 

-0.178558 

6.763176 

-0.123631 

6.772415 

-0.068047 

6.781654 

0.002376 

6.790894 

0.064782 

6.800133 

0.122427 

6.809372 

0.166800 

6.818611 

0.191209 

16.827850 

0.198699T 

6.837090 

C. 181139 

6.846329 

0.148974 

6.855568 

0.095020 

6.864807 

0.036073 

6.874046 

-0.034133 

6.883286 

-0.098084 

6.892525 

-0.157116 

6.901764 

-0.202599 

6.911003 

-0.227852 

6.920242 

-0.2372731 

6.929482 

-0.219966 

6.938721 

-0.188992 

6.947960 

-0.136223 

6.957199 

-0.076677 

6.966438 

-0.007875 

6.975677 

0.056634 

6.984917 

0.116522 

6.994156 

0.160240 

7.003395 

0.189328 

7.012634 

0.1947531 

7.021873 

0.182232 

7.031113 

0.147074 

7.040352 

0.097537 

7.049591 

0.035587 

7.058830 

-0.032350 

7.068069 

-0.096843 

7.077309 

-0.157209 

7.086548 

-0.199349 

7.095787 

-0.228534 

17. 105026 

-0.2325351 

7.114265 

-0.218630 

7.123505 

-0.183318 

7.132744 

-0.131329 

7.141983 

-0.069772 

7.151222 

-0.000209 

7.160461 

0.065476 

7.169701 

0.124636 

7.178940 

0.171080 

7.188179 

0.195974 

7.197418 

0.206602 

7.206657 

0.187360 

7.215897 

0.158368 

7.225136 

0.102800 

7.234375 

0.044973 

7.243614 

-0.025200 

7.252853 

-0.089960 

7.262093 

-0.148125 

7.271332 

-0.194823 

7.280571 

-0.219536 

7.289810 

-0.229761 

7.299049 

-0.212095 

7.308289 

-0.181532 

7.317528 

-0.127805 

7.326767 

-0.068328 

Tabic  4.1  Data  for  backscartercd  wave  under  coincident  condition  with  xo  =  34.0. 


59 


T 

P 

6.365273 

-0.023188 

6.381806 

-0.439772 

6.398339 

-0.706256 

6.414872 

-0.700890 

6.431405 

-0.446606 

6.447938 

-0.021616 

6.464471 

0.399340 

6.481004 

0.670664 

6.497537 

0.672829 

6.514070 

0.416703 

6.530602 

-0.011839 

6.547135 

-0.438064 

6.563668 

-0.707348 

6.580201 

-0.713637 

6.596734 

-0.455063 

6.613267 

-0.034457 

6.629800 

0.393084 

6.646333 

0.658783 

6.662866 

0.665960 

6.679399 

0.410796 

6.695931 

-0.012120 

6.712464 

-0.437461 

6.728997 

-0.705638 

6.745530 

-0.713730 

6.762063 

-0.457097 

6.778596 

-0.036101 

6.795129 

0.390068 

6.811662 

0.659242 

6.828195 

0.667313 

6.844728 

0.413183 

6.861260 

-0.005725 

6.877793 

-0.431122 

6.894326 

-0.699590 

6.910859 

-0.707751 

6.927392 

-0.453401 

6.943925 

-0.033302 

6.960458 

0.393116 

6.976991 

0.662974 

6.993524 

0.672867 

7.010056 

0.419696 

7.026589 

0.002894 

7.043122 

-0.426468 

7.059655 

-0.694408 

T 

P 

6.373540 

-0.238160 

6.390073 

-0.597940 

6.406606 

-0.741619 

6.423139 

-0.604760 

6.439672 

-0.246643 

6.456204 

0.203273 

6.472737 

0.564193 

6.489270 

0.705618 

6.505803 

0.571213 

6.522336 

0.212836 

6.538869 

-0.232127 

6.555402 

-0.600266 

6.571935 

-0.7434291 

6.588468 

-0.615201 

6.605000 

-0.257137 

6.621533 

0.188746 

6.638066 

0.553324 

6.654599 

0.6990541 

6.671132 

0.565667 

6.687665 

0.209122 

6.704198 

-0.236750 

6.720731 

-0.602479 

6.737264 

-0.7458431 

6.753797 

-0.614909 

6.770329 

-0.258647 

6.786862 

0.187598 

6.803395 

0.551739 

6.819928 

0.6970351 

6.836461 

0.568729 

6.852994 

0.214244 

6.869527 

-0.229921 

6.886060 

-0.593886 

6.902593 

-0.740077 

6.919126 

-0.611882 

6.935658 

-0.257104 

6.952191 

0.188291 

6.968724 

0.554449 

6.985257 

0.702523 

7.001790 

0.577343 

7.018323 

0.222881 

7.034856 

-0.220533 

7.051389 

-0.587425 

7.067922 

-0.734579 

Table  4.2  Data  for  backscattered  wave  under  coincident  condition  with  xq  55  38.0. 
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T 

P 

T 

P 

6 . 350278 

0.559231- 

6.357758 

1.161832 

6.365237 

1.646473 

6.372717 

1.967741 

6.380197 

2.097844 

6.387676 

2.019434 

6.395156 

1.732595 

6.402636 

1.275286 

6.410115 

0.689955 

6.417595 

0.042314 

6.425075 

-0.617974 

6.432554 

-1.220983 

6.440034 

-1.70833/ 

6.447514 

-2.024584 

6.454993 

-2.14501* 

6.462473 

-2.060061 

6.469953 

-1.777786 

6.477432 

-1.320424 

6.484912 

-0.730696 

6.492392 

-0.079298 

6.499871 

0.582484 

6.507351 

1.179777 

6.514831 

1.670031 

6.522310 

1.985994 

6.529790 

2.108430 

6.537270 

2.020617 

6.544749 

1.737244 

6.552229 

1.280616 

6.559709 

0.691401 

6.567188 

0.039001 

6.574668 

-0.622499 

6.582148 

-1.221521 

6.589627 

-1.709781 

6.597107 

-2.025415 

6.604587 

-2.147242 

6.612066 

-2.060807 

6.619546 

-1.774627 

6.627026 

-1.316886 

6.634505 

-0.727802 

6.641985 

-0.075336 

6.649465 

0.585943 

6.656944 

1.189285 

6.664424 

1.671964 

6.671904 

1.991014 

16.679383 

2.1095T9T 

6.686863 

2.026320 

6.694343 

1.737120 

6.701822 

1.275133 

6.709302 

0.690124 

6.716782 

0.034059 

6.724261 

-0.627559 

6.731741 

-1.234824 

.(L73.9221 

-1.717731 

6.746700 

-2.035459 

16.754180 

-2.1584581 

6.761660 

-2.069981 

6.769139 

-1.785295 

6.776619 

-1.325667 

6.784099 

-0.735967 

6.791578 

-0.078363 

6.799058 

0.587846 

6.806538 

1.188315 

6.814017 

1.675397 

6.821497 

1.995219 

6.828977 

2. 1168501 

6.836456 

2.034254 

6.843936 

1.746574 

6.851416 

1.286531 

6.858895 

' 0.699010 

6.866375 

0.039062 

6.873855 

-0.623633 

6.881334 

-1.229516 

6.888814 

-1.714690 

6.896294 

-2.033982 

16.903773 

-2. 1602281 

6.911253 

-2.072351 

6.918733 

-1.788444 

6.926212 

-1.330139 

6.933692 

-0.741370 

6.941172 

-0.084951 

6.948651 

0.580860 

6.956131 

1.187328 

6.963611 

1.674013 

6.971090 

1.997989 

6.978570 

2.120455 

6.986050 

2.036927 

6.993529 

1.750979 

7.001009 

1.290996 

Table  4.3  Data  for  backscattered  wave  under  coincident  condition  with  xq  *  42.0. 


T 

P 

T 

P 

6.597102 

-0.025907 

6.603931 

-0.973626 

6.610761 

-1.825780 

6.617590 

-2.503560 

6.624419 

-2.936231 

6.631248 

-3.079511 

6.638078 

-2.929365 

6.644907 

-2.493511 

6.651736 

-1.812736 

6.658566 

-0.959201 

6.665395 

-0.012516 

6.672224 

0.933062 

6.679053 

1.785619 

6.685883 

2.464755 

6.692712 

2.897599 

6.699541 

3.043399 

6.706370 

2.893751 

6.713200 

2.455397 

6.720029 

1.773091 

6.726858 

0.918041 

6.733687 

-0.031260 

6.740517 

-0.978651 

6.747346 

-1.832103 

6.754175 

-2.510407 

6.761004 

-2.942909 

'6.767834 

-3.0896631 

6.774663 

-2.939170 

6.781492 

-2.501066 

6.788321 

-1.819549 

6.795151 

-0.964964 

6.801980 

-0.016426 

6.808809 

0.932137 

6.815639 

1.786369 

6.822468 

2.464557 

6.829297 

2.899034 

6.836126 

3.0478141 

6.842956 

2.897290 

6.849785 

2.459860 

6.856614 

1.779217 

6.863443 

0.924893 

6.870273 

-0.024116 

6.877102 

-0.973685 

6.883931 

-1.826900 

6.890760 

-2.505418 

6.897590 

-2.941862 

16.904419 

-3 . 089979f 

6.911248 

-2.938947 

6.918077 

-2.503733 

6.924907 

-1.823455 

6.931736 

-0.967920 

6.938565 

-0.020529 

6.945395 

0.928623 

6.952224 

1.783107 

6.959053 

2.461402 

6.965882 

2.898016 

16.972712 

3 .0470171 

6.979541 

2.896656 

6.986370 

2.462214 

6.993199 

1.782508 

7.000029 

0.927581 

7.006058 

-0.018170 

7.013687 

-0.966861 

7.020516 

-1.820955 

7.027346 

-2.497132 

7.034175 

-2.933331 

7.041004 

-3.032415 

7.047833 

-2.930409 

7.054663 

-2.495530 

7.061492 

-1.816685 

7.068321 

-0.960380 

7.075150 

-0.013899 

7.081980 

0.934126 

7.088809 

1.790331 

7.095638 

2.468091 

7.102468 

2.904951 

7.109297 

3.057021 

7.116126 

2.907706 

7.122955 

2.474592 

7.129785 

1.799615 

7.136614 

0.947236 

Table  4.4  Data  for  backscattered  wave  under  coincident  condition  with  xq  =  46.0. 


T 

P 

T 

P 

6.504855 

0.011539 

6.509044 

-0.436678 

6.513232 

-0.867708 

6.517421 

-1.262138 

6.521609 

-1.602292 

6.525798 

-1.873637 

6.529986 

-2.064192 

6.534175 

-2.164548 

6.538363 

-2.169625 

6.542552 

-2.080113 

6.546741 

-1.901734 

6.550929 

-1.643118 

6.555118 

-1.314392 

6.559306 

-0.927864 

6.563495 

-0.499779 

6.567683 

-0.050313 

6.571872 

0.399248 

6.576060 

0.829382 

6.580249 

1.222500 

6.584437 

1.561802 

6.588626 

1.831941 

6.592814 

2.021164 

6.597003 

2.121951 

losnsn 

2.129898T 

6.605380 

2.043390 

6.609569 

1.865207 

6.613757 

1.603903 

6.617946 

1.272740 

6.622134 

0.886852 

6.626323 

0.461651 

6.630511 

0.013940 

6.634700 

-0.436558 

6.638888 

-0.868439 

6.643077 

-1.262070 

6.647265 

-1.601197 

6.651454 

-1.871750 

6.655643 

-2.061283 

6.659831 

-2.161029 

16. 664020 

-2.1675461 

6.668208 

-2.081394 

6.672397 

-1.905348 

6.676585 

-1.645490 

6.680774 

-1.313410 

6.684962 

-0.925579 

6.689151 

-0.500037 

6.693339 

-0.054200 

6.697528 

0.394307 

6.701716 

0.826047 

6.705905 

1.220548 

6.710093 

1.559359 

6.714282 

1.828331 

6.718471 

2.017086 

6.722659 

2.117488 

re  .726848 

2.124104] 

6.731036 

2.036243 

6.735225 

1.858479 

6.739413 

1.598993 

6.743602 

1.268452 

6.747790 

0.880678 

6.751979 

0.453133 

6.756167 

0.005424 

6.760356 

-0.442877 

6.764544 

-0.873113 

6.768733 

-1.267176 

6.772922 

-1.607440 

6.777110 

-1.877936 

6.781299 

-2.066348 

6.785487 

-2.165199 

£.789676 

-2. 17124  8[ 

6.793864 

-2.084203 

6.798053 

-1.906754 

6.802241 

-1.645953 

6.806430 

-1.313712 

6.810518 

-0.925399 

6.814807 

-0.498067 

6.818995 

-0.049818 

6.823184 

0.400200 

6.827373 

0.832200 

Table  4.5  Data  for  backscattered  wave  under  coincident  condition  with  Xo  =  50.0. 
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T 

P 

T 

P 

6.382811 

0.035442 

6.386688 

0.209292 

6.390566 

0.372359 

6.394444 

0.518352 

6.398321 

0.641011  • 

6.402199 

0.735356 

6.406076 

0.798002 

6.409954 

0.826164 

6.413832 

0.817231 

6.417709 

0.769973 

6.421587 

0.686021 

6.425465 

0.569997 

6.429342 

0.428413 

6.433220 

0.268459 

6.437098 

0.097275 

6.440975 

-0.078245 

6.444853 

-0.251158 

6.448730 

-0.414271 

6.452608 

-0.560668 

6.456486 

-0.684435 

6.460363 

-0.780417 

6.464241 

-0.843409 

,6.468119 

-0.868711] 

6.471996 

-0.854340 

6.475874 

-0.802352 

6.479752 

-0.717329 

6.483629 

-0.603766 

6.487507 

-0.465459 

6.491385 

-0.307072 

6.495262 

-0.135321 

6.499140 

0.041802 

6.503017 

0.216520 

6.506895 

0.381673 

6.510773 

0.530038 

6.514650 

0.654073 

6.518528 

0.747090 

6.522406 

0.805060 

6.526283 

0.8271531 

6.530161 

0.814186 

6.534039 

0.766641 

6.537916 

0.684721 

6.541794 

0.570449 

6.545671 

0.429062 

6.549549 

0.267968 

6.553427 

0.094795 

6.557304 

-0.083125 

6.561182 

-0.258108 

6.565060 

-0.422186 

6.568937 

-0.568257 

6.572815 

-0.690728 

6.576693 

-0.784984 

6.580570 

-0.846725 

53.4.4-4 

-0.8722681 

6.588326 

-0.859543 

6.592203 

-0.808948 

6.596081 

-0. 723501 

6.599958 

-0.608078 

6.603836 

-0.468155 

6.607714 

-0.309112 

6.611591 

-0.136769 

6.615469 

0.041721 

6.619347 

0.21E022 

6.623224 

0.383811 

6.627102 

0.531717 

6.630980 

0.655401 

,6.634857 

0.749476. 

6.638735 

0.809840  ] 

6.642612 

0.834087L 

6.646490 

0.821553 

6.650368 

0 . 773008 

6.654245 

0.690296 

6.658123 

0.576243 

6.662001 

0.435166 

6.665878 

0.273601 

6.669756 

0.100009 

6.673634 

-0.077046 

6.677511 

-0.250588 

6.681389 

-0.414793 

6.685266 

-0.563140 

6.689144 

-0.687918 

6.693022 

-0.782125 

6.696899 

-0.841541 

6.700777 

-0.864707 

6.704655 

-0.851418 

6.708532 

-0.801931 

6.712410 

-0.717605 

6.716288 

-0.601830 

6.720165 

-0.460004 

6.724043 

-0.298634 

6.727921 

-0.124493 

Table  4.6  Data  for  backscattered  wave  under  coincident  condition  with  xq  =  54.0. 


T 

P 

6.074376 

-0.004603 

6.083041 

-0.138927 

6.091706 

-0.247784 

6.100371 

-0.297811 

6.109036 

-0.281832 

6.117702 

-0.197849 

6.126367 

-0.070412 

6.135032 

0.071267 

6.143697 

0.192945 

6.152362 

0.261705 

6.161027 

0.266137 

6.169692 

0.199254 

6.178357 

0.083878 

6.187022 

-0.056422 

6.195687 

-0.184968 

6.204352 

-0.271104 

16.213017 

-0 . 294803J 

6.221683 

-0.245890 

6.230348 

-0.138985 

6.239013 

0.004546 

6.247678 

0.140881 

6.256343 

0.237589 

16.265008 

0.2627241 

6.273673 

0. 21.7343 

6.282338 

0.115705 

T 

P 

6.078709 

-0.073365 

6.087374 

-0.198649 

6.096039 

-0.281372 

16. 104704 

-0.298039) 

6.113369 

-0.247893 

6.122034 

-0.137058 

6.130699 

-0.000069 

6.139364 

0.137669 

6.148029 

0.234435 

16.156694 

0.2731101 

6.165359 

0.240549 

6.174025 

0.146085 

6.182690 

0.015065 

6.191355 

-0.124831 

6.200020 

-0.234247 

6.208685 

-0.292394 

6.217350 

-0.278369 

6.226015 

-0.199158 

6.234680 

-0.068894 

6.243345 

0.075609 

6.252010 

0.196614 

6.260675 

0.259815 

6.269341 

0.247984 

6.278006 

0.172367 

6.286671 

0.051692 

Tabic  4.7  Data  for  backscattercd  wave  under  coincident  condition  with  xq  *  58.0. 
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T 

P 

5.962266 

-0.015721 

5.970372 

-0.077883 

5.978478 

-0.122413 

5.986585 

-0.143257 

5.994691 

-0.139266 

6.002797 

-0.106693 

6.010903 

-0.060525 

6.019010 

-0.003971 

6.027116 

0.045946 

6.035222 

0.080809 

16.043328 

0.093645| 

6.051435 

0.079337 

6.059541 

0.045096 

6.067647 

-0.005789 

6.075753 

-0.057588 

6.083859 

-0.100811 

6.091966 

-0.125670 

6.100072 

-0.122697 

6.108178 

-0.098435 

6.116284 

-0.051724 

6.124391 

0.001606 

6.132497 

0.050305 

6.140603 

0.085678 

6.148709 

0 . 093451] 

6.156816 

0.078462 

6.164922 

0.039610 

T 

P 

5.966319 

-0.047572 

5.974425 

-0.103375 

5.982532 

-0.135460 

15.990638 

-6.145075] 

5.998744 

-0.125686 

6.006850 

-0.084866 

6.014956 

-0.033032 

6.023063 

0.023066 

6.031169 

0.065132 

6.039275 

0.090985 

6.047381 

0.089184 

6.055488 

0.064762 

6.063594 

0.020878 

6.071700 

-0.032454 

6.079806 

-0.080598 

6.087913 

-0 . 11655* 

16.096019 

-0. 1273l'9i 

6.104125 

-0.113186 

<5.112231 

-0.077523 

6.120337 

-0.024494 

6.128444 

0.026482 

6.136550 

0.071073 

6.144656 

0.092797 

6.152762 

0.088778 

6.160869 

0.061804 

6.168975 

0.014311 

Tabic  4.8  Data  for  backscattered  wave  under  coincident  condition  with  xo  =  62.0. 
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X 

Amplitude  by 
ray  model 

Amplitude  by 
modified  ray  model 

Amplitude  by 
PWS 

34.0 

0.169133 

0.140835 

0.2158 

38.0 

0.768413 

0.667904 

0.7213 

42.0 

2.132889 

1.919801 

2.1363 

46.0 

3.167944 

2.930751 

3.0686 

50.0 

2.394250 

2.260151 

2.1482 

54.0 

1.088952 

1.041594 

0.8506 

58.0 

0.433748 

0.418051 

0.2822 

62.0 

0.195738 

0.189482 

0.1149 

Table  5.  The  amplitudes  of  a©  wave  for  SS304  shell  by  ray  methods  and  PWS  calculation. 


Fig.  12a.  Comparison  of  the  amplitude  for  the  ray  model  and  PWS  calculation.  The  solid 
curve  is  for  the  ray  model  calculation  which  is  based  on  the  Eq.(5.5);  while  the  dots  are  for 
the  PWS  computation. 


Fig.  12b.  Comparison  of  the  amplitude  for  the  ray  model  and  PWS  calculation.  The  solid 
curve  is  for  the  ray  model  calculation  which  is  based  on  a  modification  to  Eq.(5.5)  in 
which  c/c /  is  omitted;  while  the  dots  are  for  the  PWS  computation.  Comparison  of  Fig.  12a 
and  12b  shows  that  the  exact  ray  model  is  somewhere  between  these  two. 


computations.  The  large  second  echo  visible  in  Fig.  1  Id  is  seen  to  be  the  contribution  from 
the  a©'  ray  which  has  circumnavigated  only  the  backside  of  the  shell.  Excellent  agreement 
of  the  results  in  terms  of  this  simple  calculation  is  a  manifestation  of  the  validity  of  the  ray 
model.  Of  course,  So  can  still  be  seen  especially  in  the  case  of  short  bursts  where  we  have 
more  low  frequency  components  as  presented  in  the  appendix  A.  In  addition,  the  overlap  of 
weak  ao  wave  on  ao'  can  also  be  seen  at  high  frequencies. 
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shells 

Next,  we  turn  our  attention  to  the  resonance  backscattering  of  shells.  Consider  the 
longitudinal  resonance  conditions1  ^ 


kLh  =  nrr 

(n  s  1,2,...) 

(5.6) 

kLh  =  (n+l/2)rt 

(n  *  0,1,...)  , 

(5.7) 

where  h  is  the  thickness  of  the  shell.  These  conditions  are  obtained  from  the  requirement  of 
constructive  interference  of  consecutive  internally  reflected  waves  upon  transmission  back 
into  the  water.  These  reflections  are  shown  in  Fig.7b.  Eq.(5.6)  corresponds  to  the  usual 
case  considered  here  where  the  acoustic  impedance  of  the  elastic  material  pec/  is  greater  than 
the  impedance  of  the  water  pc  while  Eq.(5.7)  is  valid  for  the  case  where  p^C/  is  less  than 
pc.  First  assume  an  ideal  fluid  shell  whose  physical  constants  are  the  same  as  the  real 
aluminum  shell  (to  be  called  aluminum/fluid  shell)  and  are  listed  in  Table  6b.  By 
investigating  such  an  ideal  object,  we  are  able  to  get  rid  of  the  effects  due  to  elastic  surface 
waves  and  transverse  waves  and  concentrate  on  the  longitudinal  resonance  echoes.  The 


Material 

p  (g/cm3) 

ci  (km/s) 

cs(km/s) 

Aluminum 

2.70 

6.42 

3.04 

Water 

1.00 

1.4825 

Tabic  6a  Material  parameters  for  Aluminum  shell  and  water. 


Fluid/water 

pe  (g/cm3) 

c  /  (km/s) 

c  (km/s) 

p  (g/cm3) 

Aluminum/water 

2.70 

6.42 

1.4825 

1.00 

SS440C/water 

7.84 

5.854 

1.479 

1.00 

Table  6b.  Material  parameters  for  fluid  shells. 
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form  functions  of  a  4%-thick  aluminum/fluid  shell  is  shown  in  Fig.  13  for  c/  ■  6.42  km/s 
where  the  form  function  for  a  real  aluminum  shell  is  also  plotted  for  comparison.  Eq.(5.6) 
predicts  a  longitudinal  resonance  at  XLr  *  340.  The  minimum  in  both  cases  is  a 
manifestation  of  this  longitudinal  resonance.  This  dip  has  been  modeled  in  detail  by  a  ray 
method11.  The  cause  of  the  large  peak  near  x  =  302  for  the  real  aluminum  shell  has  not 
been  well  understood  yet.  It  has  been  explained  as  the  prompt  radiation  in  Ref.20  and 
described  as  a  "thickness  quasiresonance"  in  Ref.17.  The  form  function  for  a  16.2%-thick 
SS440c/fluid  shell  is  shown  in  Fig.14.  This  time  Eq.(5.6)  predicts  a  longitudinal  resonance 
at  xlt-  76.8.  Again  a  dip  near  x  -  76.8  in  Fig.14  corresponds  to  the  presence  of  the 

longitudinal  resonance.  The  computed  backscattered  waves  for  the  aluminum/fluid  shell  at 
xo  *  302  and  340  are  both  shown  in  Fig.  15.  Various  burst  lengths  have  been  tested. 
Presented  in  Fig.  15a  and  15b  are  the  case  of  20-cycle  burst.  For  reference,  some 
additional  figures  are  presented  in  the  Appendix  B.  The  results  are  just  expected.  For  xq  = 
340,  the  decay  of  the  echoes  immediately  following  the  specular  reflection  is  the 
manifestation  of  a  longitudinal  resonance.  The  form  function  predicts  a  minimum  amplitude 
of  0.17,  as  incident  burst  gets  longer,  the  backscattered  echo  does  reach  the  limit  as  shown 
in  Fig.  15c  where  a  40-cycle  burst  has  been  used.  The  measured  amplitude  is  0.17  which 
agree  with  the  form  function  prediction  and  the  region  of  measurement  has  been  marked. 
We  don't  see  anything  anomalous  from  the  xQ  =  302  plot,  which  demonstrates  that  the 
anomaly  around  x  *  302  in  the  elastic  form  function  is  not  the  direct  effect  of  the 
longitudinal  resonance.  We  then  computed  the  backs cattering  for  the  SS440c/fluid  shell  for 
x0  =  71  and  76.8.  The  reason  we  are  interested  in  the  case  of  xq  =  71  will  become  clear 
later  in  Section  5.5.  The  results  are  Fig.  16  and  again  are  fully  expected.  The  decay  of  the 
echoes  corresponds  to  the  longitudinal  resonance  which  is  not  responsible  for  any  possible 
anomalous  behavior  around  x  =  71.  As  in  the  case  of  Aluminum/fluid  shell,  we  checked  the 


Fig.  15c.  Computational  backscattering  from  a  4 %- thick  aluminum/fluid  shell  for  a 
40-cycle  sine  burst  with  xq  *  340  The  measured  amplitude  as  marked  agrees  with 
the  form  function  prediction. 
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minimum  amplitude  which  is  predicted  by  the  form  function  to  be  0.63.  The  near  steady 
state  amplitude  measured  from  Fig.  16b  is  0.64  which  is  very  close  to  the  form  function 
prediction. 


5.5  Reverberations  of  longitudinal  waves  and  resonance  scattering-elastic 
shells 


Next  we  investigate  the  on  resonance  and  below  resonance  backscattering  for  a  real 
aluminum  shell.  Results  for  a  20-cycle  incident  burst  for  xo  »  302  and  340  are  depicted 
in  Fig.  17.  The  buildup  and  decay  at  the  very  early  stage  of  the  backscattered  wave  under 
the  constructive  and  destructive  resonance  conditions  are  both  evident  The  form  function 
predicts  a  maximum  amplitude  of  4.26  times  that  of  the  incident  bunts.  We  do  see  that  as 
input  burst  gets  longer,  the  amplitude  of  the  backscattering  reaches  the  limit  as  shown  in 
Fig.  17c. 

The  form  function  calculated  for  the  SS440c  shell  using  the  current  ray  model  fails 
to  recover  the  form  function  from  the  exact  PWS  calculation  in  the  x  range  of  70-74  20 . 
One  possible  cause  of  the  anomalous  behavior  near  x  =  71  is  due  to  the  fact  that  the  phase 
velocity  becomes  large  for  the  sj  Lamb  wave,  hence  the  trace- velocity  matching  angle  0/ 
defined  in  Eq.(5.1)  becomes  small  and  the  group  velocity  for  the  si  Lamb  wave  becomes 
negative  around  x  =*  71.  It  has  been  suggested  (P.L.Marston,  private  communication)  that 
the  energy  of  the  burst  is  taking  a  short  path  and  travelling  along  the  path  of  BC  in  Fig.7a 
directly  without  circumnavigating  the  shell  at  all.  This  is  one  ray  mechanism  for  producing 
prompt  radiation,  which  was  not  included  in  the  analysis  in  Ref.20.  The  computed 
backscattered  waves  from  a  real  SS440C  shell  for  xo  =  71.0  and  76.8  are  shown  in  Fig.  18 
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which  is  far  more  complicated  than  the  case  of  a  SS440c/fluid  shell.  Difference  between 
elasdc  and  fluid  shell  results  is  rapidly  evident  in  the  time  domain.  It  is  seen  by  x-  1 . 5 , 

which  means  there  is  a  prompt  mechanism  for  scattering  in  the  elastic  shell  not  included  in 
the  fluid  case.  We  also  see  that  the  prompt  radiation  may  be  out  of  phase  with  the  specular 
reflection  therefore  the  effect  can  be  either  constructive  or  destructive,  which  is  the  case  of 
aluminum  shell  and  SS440C  shell  respectively. 
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Chapter  6 

Discussion  and  Conclusion 

The  canonical  problem  of  the  backscattering  of  a  short  tone  burst  from  elastic 
spherical  shells  has  been  computationally  modeled  and  transient  scattered  waves  has  been 
calculated.  The  algorithms  used  are  simple  trapezoidal  method  to  calculate  the  Fourier 
integral  and  FFT  routines.  Both  methods  have  turned  out  to  be  equally  successful,  except 
that  the  physical  meanings  are  more  clear  in  the  former  case  while  the  latter  can  do  the 
transformations  much  faster.  Most  results  presented  in  the  report  are  from  the  calculation  of 
FFT  routines.  No  fundamental  complications  is  in  this  approach.  Any  imperfection  in  the 
calculation  appears  to  be  due  to  the  necessity  of  assuming  a  finite  bandwidth.  This  was 
taken  to  be  much  larger  than  the  carrier  frequency  for  all  cases  shown. 

Previous  measurements  by  Williams  and  Marston  for  the  tungsten  carbide  sphere 
and  by  S.  Kargl  and  P.L.Marston  for  the  16.2%-thick  elastic  SS440c  shell  have  been 
confirmed.  Along  with  Fourier  Transform  calculations,  a  ray  model  of  Sun  and  Marston19 
for  the  backscattering  from  thin  shells  has  been  carefully  investigated.  The  results 
demonstrate  that  the  thin  elastic  shells  like  5%-thick  SS440c  and  2.5%-thick  SS304  shells 
support  ao'  wave  which  can  be  "creeping"  and  "trapped"  waves.  The  phase  and  group 
velocities  for  ao'  are  both  close  to  the  sound  speed  in  the  surrounding  water.  This  particular 
surface  wave  contributes  to  the  form  function  in  the  form  of  a  big  hump  or  mid-frequency 
enhancement  for  thin  shells.  The  central  frequency  of  the  hump  varies  inversely  with  the 
thickness  of  the  shell  but  the  backscattered  echoes  due  to  the  effect  have  an  almost  constant 
arrival  time  (in  units  of  a/c)  relative  to  the  specular  reflection.  In  the  system  point  of  view. 


the  scatterer  is  like  a  band-pass  filter.  The  frequency  components  in  the  range  of  the  hump 
are  relatively  enhanced.  Because  of  its  large  magnitude,  this  enhancement  and  other 
properties  of  the  associated  surface  wave  may  be  of  great  practical  use.  Appendix  C  shows 
that  the  enhancement  also  contributes  in  a  distinct  way  to  the  sphere’s  impulse  response. 

The  thickness  resonance  conditions  are  also  examined  by  computing  the 
backscattering  under  those  conditions.  For  the  fluid  stainless  steel  and  aluminum  shell 
cases,  the  agreement  between  the  calculated  backscattering  and  the  prediction  by  PWS 
form  function  under  the  longitudinal  thickness  resonance  is  good.  Since  the  fluid  shell  case 
has  been  well  modeled  by  a  ray  method,  this  also  indicates  a  good  agreement  with  the  ray 
prediction.  A  prompt  radiation  that  is  currently  being  modeled  appears  to  be  important  in 
elastic  shell  cases  and  the  amplitude  of  the  computed  backscattering  for  the  aluminum  shell 
does  reach  the  limit  predicted  by  the  PWS  form  function.  The  case  of  a  thick  elastic 
SS440C  shell  (16.2%-thick)  is  more  complicated  in  that  the  prompt  radiation  appears  to 
interfere  destructively  with  the  specular  reflection. 


oooo 


\ 


87 


Appendix 


A.  List  of  programs 


Al.  Program  computing  the  form  function  of  elastic  spheres 


ccccccccccccccccccccccccccccccccccccccccccccccccccccccccc 

c 

THIS  PROGRAM  CALCULATES  FORM  FUNCTIONS  USING  CHIVERS  C 
PROGRAMS .  C 

C 

CCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCCC 

c 

IMPLICIT  DOUBLE  PRECISION  (P) 

DIMENSION  PB{800) ,PCB(800) ,PSB(800) 

C  REAL  RAM 

COMPLEX  CFAC 

DATA  RHO, VC, VS, VW/13. 80, 6860. 0,4 185. 0,1476.0/ 
c  data  rho,vc,vs,vw/7.70,  5960.0,3240.0,1410.0/ 

RH=1.0/(RHO*2.0) 

RAT1=  VW/VC  , 

RAT2=  VW/VS 

XMAX=*  400.0 
X0=0.0 
K=0 

10  k=k+l 

x=x0+k*0 . 05 
T=1 . 25*X 
X1=RAT1*X 
X2=RAT2*X 
X2S=X2**2 
XSN=SIN (X) 

XCS=COS(X) 

CALL  BESS (X, PB,T) 

TERM=  (XSN/X) /PB ( 1) 

BN1=  PB ( 2 ) *TERM 
CALL  BESS (XI , PCB, T) 

CBN1=  PCB ( 2) /PCB ( 1) 

CALL  BESS (X2 , PSB, T) 

SNN0=-XCS/X 
XSIGN=1 . 0 
SBET=0 . 0 
SALP=0 . 0 
N=0 

SUM=0 . 0 

CALL  STVAL (X , XI , X2 , ALPN , BETN , CBN1 , BN1 , SNN1 , RH , XSN , X2S , XCS ) 
20  CONTINUE 

CBN 2=  PCB (N+3 ) /PCB (N+2 ) 

Q=FLOAT (N+l) 


nno 
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ANUM«(Q*CBN1) -(X1*CBN2) 

ADEN-( (Q-1.0) *CBN1) -(X1*CBN2) 

DNUM-( ( (X2S/2.0)-Q* (Q-1.0) ) *CBN1) - (2 . 0*X1*CBN2 ) 
SBN1-1.0 

SBN2-  PSB (W+3 ) /PSB (N+2 ) 

BNUM* ( 2 . 0*Q* (Q+l • 0 ) ) *SBN1 

EDEN“ ( ( ( 2 . 0* (Q**2) ) -(X2S+2.0) ) *SBN1) + (2 . 0*SBN2*X2) 
ENUM=*2 .0*Q*(Q+1.0)  *  {  (  (1.0-Q)  *SBN1) + (X2  *SBN2 )  ) 
FN1=X2S*RH*  (  ( ANUM/ADEN)  -  ( BNUM/EDEN)  )  /  ( (DNUM/ADEN) 

&- (ENUM/EDEN) ) 

30  CONTINUE 

BN2»  PB(N+3)*TERM 
GN1= ( (FN1-Q) *BN1) + (X*BN2 ) 

BN1=BN2 

SNN2* ( ( (2 . 0*Q) +1.0) *SNN1/X) -SNNO 
HN1»( (FN1-Q) *SNN1) + (X*SNN2 ) 

SNN0=SNN1 

SNN1-SNN2 

DEN= (GN1**2 ) + (HN1**2 ) 

ALPN1=-(GN1**2)/DEN 
BETN1=»-  (GN1*HN1 )  /DEN 

ADD=Q* (ALPN+ALPN1+ ( 2 . 0*ALPN*ALPN1)  +  ( 2 . 0  *  BETN  *  BETN 1 ) ) 
FAC=XSIGN* ( 2 . 0*Q-1 .0) 

SALP=*SALP+FAC*ALPN 
SBET=S  BET+FAC*  BETN 
ALPN=ALPN1 
BETN=B£TN1 
SUM=SUM+ADD 

IF  ( (ABS (ADD/SUM) ) .LT.1.0E-8)  GO  TO  90 

J=0 

N=N+1 

XSIGN=-XSIGN 
GO  TO  20 
90  J=J+1 

IF  (J.EQ.5)  GO  TO  100 
N=N+ 1 

XSIGN=-XSIGN 
GO  TO  20 

100  YP=(-4 . 0/ (X**2) ) *SUM 
FAC=-XSIGN* (2 . 0*Q+1 . 0) 

SALP=SALP+FAC*ALPN1 
SBET^S BET+FAC* BETN 1 
CFAC=CMPLX (SBET, -SALP) * (2 . 0/X) 

FINF=CABS ( CFAC) 

KAM (K) =X 
FABS ( K) =FINF 

WRITE (7,40)  KAM (K) , FABS (K) 

WRITE (7,40)  X , CFAC 

40  FORMAT (  3fl2.6) 

IF(X.Gt.XMAX)  GO  TO  50 
GO  TO  10 
50  CONTINUE 


OO  o  o 
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STOP 

END 


SUBROUTINE  STVAL ( X , XI , X2 , ALPN , BETN , CBN1 , BN1, SNN1 , RH , XSN , X2S 
&,XCS) 

TERM=CBN1*X1 

FN=-X2S*RH*TERM/ ( (X2S/2 . 0) -2 . 0*TERM) 

SNN1=- ( (XCS/ (X**2 ) ) + (XSN/X) ) 

GN=(FN* (XSN/X) )+(BNl*X) 

HN=(-(FN*XCS)/X)+(X*SNN1) 

DEN= (GN**2)+(HN**2) 

ALPN=-(GN**2)/DEN 
BETN=- (GN*HN) /DEN 
RETURN 
END 


SUBROUTINE  BESS (X, PA, T) 

IMPLICIT  DOUBLE  PRECISION  (P) 

DIMENSION  PA (8 00) 

L=IFIX (T) +15 
PA (L+2 ) =0 . 0 
PA (L+l) =1 . 0E-30 
DO  10  1=1 , L 
M= (L+l) -I 
R=FLOAT (M) 

PA (M) = ( ( ( 2 . 0*R) +1.0) *PA(M+l)/X)-PA(M+2) 
10  CONTINUE 
RETURN 
END 


90 


A2.  Program  computing  the  form  function  of  elastic  spherical  shells 


REAL+8  X, XT, XL, QT, QL, BA, RHO, RHOE , XMIN, XMAX, DX,  JX (900) , 

+  JXT(900) , JXL(900) ,JQT(900) , JQL(900) , YX (900)  ,YXT(900) , 

+  YXL( 900) , YQT ( 900) ,YQL(900) ,0(5,5) ,01(5,5)  , 

+  C , DETOR , DETOI , 

+  CT , CL, XTS , QTS , DNR , DNI , COEF , FORM , SUMR , SUMI 
INTEGER  NMAX , I , J , MAX , K , INDX ( 5 ) 

C-1.4790D0 
CT-3. 141000 
CL*5. 6750D0 
RH0=>1 . 0D0 
RHOE=*7 . 570D0 
BA*0 . 9750D0 

XMIN-478.0D0 
XMAX»48 0.000 
X=»XMIN 
DX=*0 . 050D0 

MAX=INT ( (XMAX-XMIN) /DX) 

DO  1000  1*1,  MAX 

X=X+DX 

XT*X*C/CT 

XL*X*C/CL 

QT*XT*BA 

QL*XL*BA 

QTS=QT*QT 

XTS*XT*XT 

IF  (  (0.20  .LT.  X)  .AND.  (X  .LT.  8.0))  THEN 
NMAX*2+NINT (X+4 . 0*DEXP ( DLOG ( X) / 3 . 0 ) ) 

ELSE 

NMAX*3+NINT(X+4 . 05*DEXP ( CLOG ( X) /3 .0) ) 

ENDIF 

CALL  BESSEL (NMAX+2 , X, JX) 

CALL  BESSEL (NMAX+2 , XT , JXT) 

CALL  BESSEL(NMAX+2 , XL, JXL) 

CALL  BESSEL (NMAX+2, QT,JQT) 

CALL  BESSEL (NMAX+2 , QL, JQL) 

CALL  NEUMANN (NMAX+ 2, X,YX) 

CALL  NEUMANN (NMAX+2, XL, YXL) 

CALL  NEUMANN (NMAX+2, XT, Y XT) 

CALL  NEUMANN (NMAX+2 , QL, YQL) 

CALL  NEUMANN (NMAX+ 2, QT, YQT) 

SUMR=0 . 0 
SUMI=0 . 0 
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200 


DO  500  J«l,  NMAX 
FJ-DBLE(J-l) 

D  ( 1 , 1 )  “  ( RHO/RHOE ) *XTS*JX(J) 
D(1,2)»(2.0*FJ*(FJ-1.0)-XTS) *JXL(J) 
D(1,2)=D(1,2)+4.0*XL*JXL(J+1) 

D(1,3)”(2.0*FJ*(FJ-1.0) -XTS) *YXL( J) 

D(  1 , 3 ) =D ( 1 , 3 ) +4 . 0*XL*YXL( J+l ) 

D ( 1 , 4 ) ■ (FJ-1 . 0) *JXT ( J) -XT* JXT (J+l) 
D(1,4)«2.0*FJ*(FJ+1) *D(1, 4) 

D  ( 1 , 5)  **  (FJ-1 . 0)  *YXT  ( J)  -XT*YXT  (J+l) 
D(1,5)“2.0*FJ*(FJ+1)*D(1,5) 

D(2, 1)  s»-FJ*JX(J)-rX*JX(J  +  l) 

D  ( 2 , 2 ) =*FJ*JXL{  J) -XL*JXL( J+l ) 

D ( 2 , 3 ) =FJ*YXL( J) -XL*YXL ( J+l) 

D(2, 4) =FJ* (FJ+1 . 0) *JXT ( J) 

D(2 , 5)  =*FJ*  (FJ+1. 0)  *YXT( J) 

D(3 , 1) *0 . 0D0 

D(3,2)=2.0*( (1.0-FJ) *JXL(J)+XL*JXL(J+1) ) 

D(3,3)=2.0*( (1.0-FJ) *YXL ( J) +XL*YXL( J+l) ) 

D ( 3 , 4 ) = (XTS-2 .0*(FJ**2)+2.0) *JXT ( J) -2 . 0* XT* JXT (J+l) 
D(3,5)=(XTS-2.0*(FJ**2)+2.0) *YXT ( J) -2 . 0*XT*YXT ( J+l) 

D (4 , 1) =0 . 0D0 

D(4,2)=(2.0*FJ* (FJ-1 . 0) -QTS) *JQL(J) 

D(4,2)=D(4,2)+4. 0*QL*JQL( J+l) 

D(4,3)=(2.0*FJ*(FJ-1.0) -QTS) *YQL(J) 
D(4,3)=D(4,3)+4.0*QL*YQL(J+1) 

D ( 4 , 4 ) =2 . 0*FJ* (FJ+1 . 0) *( (FJ-1 . 0) *JQT ( J) -QT*JQT ( J+l) ) 
D(4,5)=2.0*FJ*(FJ+1.0)*( (FJ-1 . 0) *YQT (J) -QT*YQT (J+l) ) 
D ( 5 , 1 ) =0 . 0D0 

D(5,2)=2.0*( (1.0-FJ) *JQL( J) +QL*JQL( J+l) ) 

D(5,3)=2.0*( (1.0-FJ) *YQL(J) +QL*YQL( J+l) ) 

D(5, 4)=(QTS-2.0*(FJ**2)+2.0) *JQT(J) -2 . 0*QT*JQT(J+1) 

D ( 5 , 5) = (QTS-2 .0*(FJ**2)+2.0) *YQT(J) -2 . 0*QT*YQT(J+1) 

DI (1, 1)= (RHO/RHOE) *XTS*YX(J) 

DI (2 , 1) =-FJ*YX ( J) +X*YX (J+l) 

DI  ( 3 , 1 ) =0 . ODO 
DI ( 4 , 1 ) =0 . ODO 
DI ( 5 , 1 ) =0 . ODO 
DO  200  K=2 ,  5 

DI ( 1 , K) =D ( 1 , K) 

DI (2 , K) =D (2 , K) 

DI ( 3 , K) =D(3 , K) 

DI ( 4 , K) =D (4 , K) 

DI ( 5 , K) =D ( 5 , K) 

CONTINUE 

CALL  LUDCMP ( D , 5 , 5 , INDX , DETDR) 

CALL  LUDCMP ( DI , 5 , 5 , INDX , DETDI ) 

DO  300  K=1 ,  5 

DETDR=DETDR*  D ( K ,  K ) 

DETDI=DETDI*DI (K , K) 


300  CONTINUE 
DNR-DETDR 
DNI-DETDI 

COEF=( (-1.0) **(J+1) ) *(2.0*FJ+1.0) 
SUMR-SUMR+COEF* (DNI/ (DNR+DNI* (DNI/DNR) ) ) 
SUMI-SUMI+COEF* (DNR/ (DNR+DNI* (DNI/DNR) ) ) 
500  _  CONTINUF 

FORM" (2 . 0/X) *DSQRT (SUMR**2+SUMI**2) 
SUMR-  SUMR* 2 . 0/X 
SUMI-  SUMI*2 . 0/X 
WRITE (7 , 2000)  X, SUMR, SUMI 
1000  CONTINUE 

2000  FORMAT (3F12 . 6) 

END 


SUBROUTINE  NEUMANN (N, X, YN) 

REAL*8  X, YN (900) 

INTEGER  I , N 

IF  (  (X  .LT.  0.0)  .OR.  (N  .LT.  2)  )  THEN 
WRITE (6,*)  'FAILURE  IN  NEUMANN  FCN  ROUTINE' 
GOTO  30 
ENDIF 

YN ( 1 ) =-DCOS (X) /X 

YN  ( 2 )  =-DCOS  (X)/(X**2)-DSIN(X)/X 

DO  10  1=1, N-l 

YN (I  +  2)  =  ( (2. 0*DBLE (I-l)+3.0)/X) *YN ( 1  +  1) -YN ( I ) 
10  CONTINUE 

30  RETURN 

END 

SUBROUTINE  BESSEL(N, X, JN) 

REAL* 8  X , JN ( 900 ) , RATIO , JO , JMX , JMX1 

INTEGER  I, MAX ,N, NT 

JMX=1. 0 

JMX 1=1 . 0D-100 

MAX=N+20 

IF  (  (X  .LT.  0.0)  .OR.  (N  .GT.  9  :0)  )  THEN 
WRITE (6,*)  'FAILURE  IN  BESSEL  FCN  ROUTINE' 
GOTO  90 
ENDIF 

J0=DSIN(X)/X 
IF  (N  .GT.  INT (X) )  THEN 
DO  10  1=1, MAX-1 
NT=MAX-I 

JN (NT) =( (2 . 0*DBLE (NT) +1 . 0)/X) *JMX-JMX1 
JMX 1= JMX 
JMX=JN (NT) 

10  CONTINUE 

RATIO=J  0/ JN ( 1 ) 

DO  20  1=1, N+l 

JN ( I ) =RATIO* JN ( I ) 
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20  CONTINUE 

ELSE 

JN ( 1) =J0 

JN ( 2 ) =DSIN (X) / (X**2 ) -DCOS (X) /X 
DO  35  1*1,  N-2 

JN ( 1+2 ) * ( (2.0*DBLE(I)+1.0)/X)*JN(I+1)-JN(I) 
35  CONTINUE 

ENDIF 

90  RETURN 

END 

SUBROUTINE  LUDCMP ( A, N, NP, INDX,  D) 

IMPLICIT  REAL* 8  (A-H,  O-Z) 

PARAMETER (NMAX=100,  TINY=1 . 0D-100) 

DIMENSION  A  (NP,  NP)  ,INDX(N)  ,  W(NMAX) 

D-1.0 

DO  12  1=1,  N 

AAMAX=0 . 0 
DO  11  J=l,  N 

IF  ( DABS ( A ( I , J) )  .GT.  AAMAX)  AAMAX=DABS (A(I , J) ) 

11  CONTINUE 

IF  (AAMAX  .EQ.  0.0)  PAUSE  'SINGULAR  MATRIX' 

W ( I ) = 1 . 0/ AAMAX 

12  CONTINUE 

DO  19  J=l,  N 
IF  (J  .GT.  1)  THEN 
DO  14  1=1,  J-l 
SUM=A ( 1 , J ) 

IF  (I  .GT.  1)  THEN 
DO  13  K=1 ,  1-1 
SUM=SUM- A ( I , K) * A ( K , J ) 

13  CONTINUE 

A(I , J) =SUM 

ENDIF 

14  CONTINUE 
ENDIF 
AAMAX* 0.0 

DO  16  I=J,  N 
SUM=A (I , J) 

IF  (J  .GT.  1)  THEN 
DO  15  K=1 ,  J-l 

SUM=SUM-A(I,K) *A(K, J) 

15  CONTINUE 

A (I, J)=SUM 
ENDIF 

DUM= W ( I ) *  DABS ( SUM ) 

IF  (DUM  .GE.  AAMAX)  THEN 
IMAX=I 
AAMAX* DUM 
ENDIF 

16  CONTINUE 

IF  (J  .NE.  IMAX)  THEN 
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DO  17  K-l,  N 

DUM«A(IMAX,K) 

A(IMAX,K)«A(J,K) 

A  ( J ,  K) “DUM 

17  CONTINUE 
D=-D 

W  { IMAX)  “W  ( J) 

ENDIF 

INDX (J) “IMAX 

_  IF  (J_.NE.  N)  THEN 

IF  (A( j, J)  . EQ.  0.0)  A(J,J)»TINY 
DUM=1 . 0/A ( J , J) 

DO  18  I-J+l,  N 
A(I,J)-A(I,J) *DUM 

18  CONTINUE 
ENDIF 

19  CONTINUE 

IF  (A (N, N)  .EQ.  0.0)  A(N,N) =TINY 

RETURN 

END 


95 


A3.  Program  computing  the  form  function  of  fluid  spherical  shells 

REAL* 8  X,XL,QL,BA,RHO,RHOE,XMIN,XMAX,DX,JX(600) , 

+  JXL(600) , JQL( 600) ,YX(600) , YXL(600) ,YQL(600) ,D(3,3) ,DI(3,3) , 

+  C , DLTDR , DETDI , CL , DNR ,  DNI , COEF , FORM , SUMR , SUMI 

INTEGER  N'MAX ,  I ,  J ,  MAX ,  X ,  ZNOX  ( 3 ) 

C-1.4790D0 

RHO-l.ODO 

RHOE-7.840DO 

CL-5.854D0 

BA-0.838D0 

XMIN-0.0D0 

XMAX-400.000 

X-XMIN 

DX-0.050D0 

PI-3.1415926 

MAX-INT ( (XMAX-XMIN) /DX) 

DO  1000  I»1,MAX 

X-I*DX 

XL-X*C/CL 

QL-XL*BA 

IF  (  (0.20  .LT.  X)  .AND.  (X  .LT.  8.03)  THEN 
NMAX-2+NINT ( X+4 . 0*DEXP(DLOG(X) /3 . 0)  ) 

ELSE  * 

NMAX-3+NINT (X+4 . 05*DEXP (DLOG (X) /3 . 0) ) 

ENDIF 

CALL  BESSEL(NMAX+2,X,JX) 

CALL  BESSEL (NMAX4- 2  ,  XL,  JXL) 

CALL  BESSEL (NMAX+2 ,QL,JQL) 

CALL  NEUMANN (NMAX+ 2, X,YX) 

CALL  NEUMANN (NMAX+2, XL, YXL) 

CALL  NEUMANN (NMAX+2, QL,YQL) 


SUMR-0 . 0 
SUMI-0.0 

DO  500  J-l,  NMAX 
FJ-DBLE(J-l) 

D( 1 , 1) - (RHO/RHOE) * (C/CL) *JX(J) 
D ( 1 , 2 ) *JXL( J) 

D(l, 3)-YXL(J) 

D(2 , 1) -FJ*JX (J) /X-JX (J+l) 

D (2 , 2 ) -FJ*JXL( J) /XL-JXL( J+l ) 

D ( 2 , 3 ) »FJ*YXL( J) /XLWYXL( J+l) 


D (3 , 1) "0 . 0D0 
D  ( 3 , 2 ) *JQL( J) 

D(3,3)»YQL(J) 

DI ( 1 , 1) * (RHO/RHOE) * ( C/CL) *YX  ( J) 

DI ( 2 , 1) *FJ  *YX ( J) /X-YX ( J+l) 

DI ( 3 , 1) *0 . 0D0 
DO  200  K-2,  3 

DI(1,K)-D(1,K) 

DI(2,K)-D(2,K) 

DI ( 3 , K) *D ( 3 , K) 

200  continue 

CALL  LUDCMP ( D , 3 , 3 , INDX , DETDR) 

CALL  LUDCMP (DI, 3, 3 , INDX, DETDI) 

DO  300  K-l,  3 

DETDR«DETDR*D(K,K) 

DETDI»DETDI*DI(K,K) 

300  CONTINUE 

DNR»DETDR 
DNI-DETDI 

COEF-( (-1.0) **(J+1) ) *(2.0*FJ+1.0) 
SUMR«SUMR+COEF*(DNI/ (DNR+DNI* (DNI/DNR) ) ) 
SUMI»SUMI+COEF* (DNR/ (DNR+DNI* (DNI/DNR) ) ) 

500  CONTINUE 

SUMR»(2.0/X) *SUMR 
SUMI= ( 2 . 0/X) *SUMI 
write  (7,2000)  x,sumr,sumi 
1000  CONTINUE 

2000  FORMAT (3F12. 6) 

END 

SUBROUTINE  NEUMANN (N, X, YN) 

REAL*8  X, YN (600) 

INTEGER  I , N 

IF  (  (X  .LT.  0.0)  .OR.  (N  .LT.  2)  )  THEN 
WRITE (6,*)  ' FAILURE  IN  NEUMANN  FCN  ROUTINE' 
GOTO  30 
ENDIF 

YN(1)— DCOS(X)/X 

YN(2)— DCOS(X)/(X**2)-DSIN(X)/X 

DO  10  1*1, N-l 

YN(I+2)®( (2. 0*DBLE (1-1) +3 . 0) /X) *YN(I+1) -YN(I) 
10  CONTINUE 

30  RETURN 

END 
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SUBROUTINE  BESSEL(N, X, JN) 

REAL*8  X, JN (600) , RATIO, JO, JMX, JMX1 

INTEGER  I, MAX, N, NT 

JMX=1 . 0 

JMX1=1 . 0D-100 

MAX=N+2  0 

IF  (  (X  .LT.  0.0)  .OR.  (N  .GT.  500)  )  THEN 
WRITE (6,*)  • FAILURE  IN  BESSEL  FCN  ROUTINE ' 

GOTO  90 
ENDIF 

J0=DSIN(X)/X 
IF  (N  .GT.  INT (X) )  THEN 
DO  10  1=1, MAX-1 
NT=MAX-I 

JN (NT) = ( ( 2 . 0*DBLE (NT) +1 . 0) /X) *JMX-JMX1 
JMX 1= JMX 
JMX=JN (NT) 

10  CONTINUE 

RATIO=J  0/ JN ( 1 ) 

DO  20  1=1, N+l 

JN ( I ) =RATIO* JN ( I ) 

20  CONTINUE 

ELSE 

JN ( 1 ) =J  0 

JN ( 2 ) =DSIN (X)/(X**2) -DCOS (X) /X  ' 

DO  35  1=1,  N-2 

JN ( 1+2 ) = ( (2 . 0*DBLE (I ) +1 . 0) /X) *JN(I+1) -JN(I) 

35  CONTINUE 

ENDIF 

90  RETURN 

END 

SUBROUTINE  LUDCMP ( A, N , NP, INDX ,  D) 

IMPLICIT  REAL*8  (A-H,  O-Z) 

PARAMETER (NMAX=1 00,  TINY=1 . 0D-150) 

DIMENSION  A(NP,NP)  ,INDX(N)  ,W(NMAX) 

D=1 . 0 

DO  12  1=1,  N 

AAMAX=0 . 0 
DO  11  J=l,  N 

IF  ( DABS ( A ( I , J) )  .GT.  AAMAX)  AAMAX=DABS ( A ( I , J ) ) 

11  CONTINUE 

IF  (AAMAX  .EQ.  0.0)  PAUSE  'SINGULAR  MATRIX' 

W(I)  =1.0/  AAMAX 

12  CONTINUE 

DO  19  J=l,  N 
IF  (J  .GT.  1)  THEN 
DO  14  1=1,  J-l 
SUM=A ( I , J) 

IF  (I  .GT.  1)  THEN 
DO  13  K=1 ,  1-1 
SUM=SUM-A { I , K) *A(K, J) 
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CONTINUE 

A(I,J)*SUM 

ENDIF 
CONTINUE 
ENDIF 
AAMAX=*0 . 0 
DO  16  I*J,  N 
SUM=A(I,J) 

IF  (J  .GT.  1)  THEN 
DO  15  K*l,  J-l 

SUM=SUM-A(I,K) *A (K, J) 
CONTINUE 
A ( I , J) =SUM 
ENDIF 

DUM*  W  f I ) *  DABS ( SUM ) 

IF  (DUM  .GE.  AAMAX)  THEN 
IMAX=I 
AAMAX* DUM 
ENDIF 
CONTINUE 

IF  (J  .NE.  IMAX)  THEN 
DO  17  K-l,  N 

DUM=*A ( IMAX , K) 

A ( IMAX , K) *A ( J , K) 

A ( J , K) *DUM 
CONTINUE 
D*-D 

W  ( IMAX )  *W  ( J ) 

ENDIF 

INDX ( J) *IMAX 
IF  (J  .NE.  N)  THEN 

IF  (A { J, J)  .EQ.  0.0)  A(J, J)*TINY 
DUM*1 . 0/A(J, J) 

DO  18  I-J+l,  N 
A ( I , J) *A ( I , J) *DUM 
CONTINUE 
ENDIF 
CONTINUE 

IF  (A(N, N)  .EQ.  0.0)  A (N , N) *TINY 

RETURN 

END 


ooon  ooooo 
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A4.  Program  computing  transient  scattering:  integration 

REAL  M,A,B,X0,X,X1,X2,  TWOPI,G,REG, IMG,  ARG, REF, IMF, REAL, 
$  IMAG,  P(5000) , BO 
INTEGER  I,N, J ,N1 
WRITE ( 5 , 10) 

10  FORMAT ( ' READ  IN  N  &  Nl') 

READ  (5,20)  N,N1 
20  FORMAT (214) 

XO-49.1 

BO-4. 

TWOPI-2*3. 1415926 
M— TWOPI *  BO/XO 
XI— 0 . 05 
T1-20*M/N1 
DO  200  J-1,N1 
P(J)-00.0 
200  CONTINUE 

DO  100  I— 1,N 
T  -  -6.*M 

C  X-I*X1 

READ  (9,500  )  X, REF, IMF 

REF— REF 
IMF-+IMF 

500  FORMAT ( 3  F12 . 6 ) 

X2-X-X0 
X2-ABS (X2) 

A-X+XO 
B-X-XO 

IF  (X2.LT.1.0E-8)  GO  TO  30 
REG- ( 1-COS (M*A) ) /A- (1-COS (M*B) ) /B 
IMG-SIN (M*A) /A-SIN (M*B) /B 
IF(REG.EQ.O.O)  GO  TO  25 
ARG  -ATAN2 (IMG,  REG) 

GO  TO  28 

25  ARG  -3.1415926/2 
28  ARG  »  ARG/3. 1415926 
GO  TO  40 

30  REG- ( 1-COS (M*A) ) /A 
IMG-SIN (M*A) /A-M 
40  IMG— IMG 

G= (REG**2+IMG**2 ) **0.5 
ARG  -  ATAN2 ( IMG ,  REG) /3 . 14 15926 
WRITE  (7,50)  X, G 

REAL-REG*REF  -IMG* IMF 
IMAG»IMG*REF  +REG*IMF 
DO  300  J-l , Nl 

P(J)=(REAL*COS(X*T)+IMAG*SIN(X*T) ) *X1+P(J) 

T-T+Tl 

300  CONTINUE 
100  CONTINUE 
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T-  -6 . 0*M-T1 
DO  400  J-1,N1 
T*T+T1 

P(J)=P(J) /TWOPI 
WRITE (8,50)  T,  P(J) 
'400  CONTINUE 
50  FORMAT ( 2F12 . 6) 

STOP 

END 


t  Program  A4.  This  program  computes  the  transient  backsc  altering  from  any  scatterer 

for  a  BO-cycle  sine  burst  with  a  carrier  frequency  xq.  A  simple  trapezoidal  method 
is  used  in  doing  the  integration  in  Eq.(3.10).  The  disk  defined  as  device  9  contains 
the  data  for  the  form  function  of  any  scatterer.  The  N  and  N1  are  parameters  which 
specify,  respectively,  the  number  of  data  for  the  form  function  and  the  number  of 
points  you  want  for  each  cycle  of  the  incident  bunt. 
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A5.  Program  computing  transient  scattering:  FFT  routines 

INTEGER  I,  N,  M,  N1 

REAL  CONST,  TWOPI,  XO/B,  UL,T,T1,M0DSEQ  , MODCOEF,KA, 

$  REF, IMF,T2 , Y 

COMPLEX  C,  COEF(34280) ,H,  SEQ(34280) 

EXTERNAL  CONST,  FFTCF,  FFTCB 
COMMON  /WORKS P/  RWKSP 
REAL  RWKSP; 206000) 

CALL  IWKIN (206000) 

XO-49.1 

B=4. 

M=INT(B*25.) 

N=INT(X0*M/B/0. 05) 

IF  (N/2.GT.0)  N-N+l 
TWOPI-  2 . 0*CONST ( ' PI ' ) 

UL— TWOPI *B/X0 
T  =0. 

T2=0 . 

Tl-UL/M 

Y  -  3*FLOAT(M)/UL 
Nl=  IFIX ( Y) +1 
DO  30  1=1, N1 
SEQ(I) =0. 

MODSEQ  -REAL ( SEQ ( I ) ) 

WRITE  (6,300)  T, MODSEQ 
T-T+Tl 
30  CONTINUE 

DO  10  I-Nl+1  , N1+1+  M 
SEQ ( I ) =  SIN (X0*T2 ) 

MODSEQ  -REAL ( SEQ (I) ) 

WRITE  (6,300)  T, MODSEQ 
T2-T2+T1 
T  -T  +T1 
10  CONTINUE 

DO  20  I— N1+2+M  ,  N 
SEQ ( I ) —0 . 

MODSEQ  -REAL (SEQ (I) ) 

IF  (T.GT.5)  GO  TO  20 
WRITE  (6,300)  T, MODSEQ 
T-T+Tl 

20  CONTINUE 

CALL  FFTCF  (N,  SEQ,  COEF) 

DO  1000  1=1,  N 

KA=FLOAT(I-l) *XO*M/B/FLOAT (N) 

MODCOEF  -CABS (COEF (I) ) /FLOAT (M) 

BETA-ATAN2 ( AIMAG ( COEF ( I ) ) ,  REAL(COEF(I) ) ) /CONST (' PI M  *180 
WRITE  (7,300)  KA, MODCOEF 
WRITE  (2,300)  KA , BETA 
1000  CONTINUE 

DO  4000  1=2,8001 


102 


READ  (9,  5000)  X, REF, IMF 
5000  FORMAT  (3F12.6) 

COEF  ( I )  =COEF  ( I )  *CMPLX  ( -REF ,  -IMF) 
4000  CONTINUE 

DO  3000  I-  8001  , N 
COEF (I)=(0.,0.) 

3000  CONTINUE 

CALL  FFTCB  (N,COEF,SEQ) 

T=-UL/ FLOAT (M) 

DO  2000  1=1, N 

SEQ ( I ) *SEQ (I) *2. /FLOAT (N) 

MODSEQ  »  REAL(SEQ(I) ) 

T=  T  +T1 

IF  (T.GE.10.)  GO  TO  2500 
WRITE (8,300)  T,  MODSEQ 

300  FORMAT ( 2 F12 . 6) 

2000  CONTINUE 
2500  STOP 
END 


Program  A5.  This  program  computes  the  transient  backscattering  from  any  scattercr 
for  a  B-cycle  sine  burst  with  a  carrier  frequency  xo.  Two  FFT  routines,  namely, 
FFTCF  and  FFTCB,  in  the  CMS  main  frame  Fortran  routines  bank  are  utilized. 
The  parameters  N  and  M  specify,  respectively,  the  total  number  of  points  for  the 
input  sequence  and  the  number  of  points  per  cycle  for  the  input  burst.  The  disk 
defined  as  device  9  contains  the  data  for  the  form  function  of  any  scattercr.  The  sign 
of  these  data  may  have  to  be  changed  depending  on  the  time  convention  the  form 
function  program  is  using. 


/ 
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A6.  Illustration  about  writing  and  using  the  programs 

The  program  computing  the  form  function  of  an  elastic  sphere  was  written  by 
Olivers  and  Anson22.  Both  the  program  computing  the  form  function  of  elastic  spherical 
shells  and  that  of  fluid  shells  were  written  by  S.  Kargl.  They  are  presented  here  just  for 
user's  convenience.  These  programs  are  found  written  in  e  '1C0t  time  convention  whose 
results  can  be  used  in  present  computation  without  any  change.  To  compute  the  form 
function  for  high  ka  values,  which  is  often  the  case  for  transient  calculations,  one  needs  to 
expand  the  dimension  of  the  arrays  defined  in  the  programs.  But  this  is  not  limitless.  These 
programs  can  not  calculate  the  form  function  for  arbitrarily  high  ka  values.  This  could  be  a 
problem  sometimes.  As  in  the  case  of  resonance  scattering  by  aluminum  shell,  the  carrier 
frequency  is  as  high  as  x<j  =  340,  but  the  form  function  can  only  be  calculated  up  to  x  = 
480  for  the  material  constants  we  are  using.  The  cut-off  high  frequency  of  480  gives  us  a 
real  incident  burst  as  shown  in  Fig.  A1  which  is  far  from  perfect.  When  interpreting  the 
transient  scattering  in  this  case,  one  has  to  take  into  the  consideration  the  no  longer  small 
difference  between  the  ideal  and  real  incident  bursts. 

The  programs  computing  the  transient  scattering  are  also  written  in  e-itot  time 
convention  but  the  FFT  routines  called  from  IMSL  bank  are  in  e+irnt  time  convention. 
Since  our  input  is  a  real  sequence,  the  difference  in  time  convention  gives  a  different  sign 
for  the  imaginary  part  of  the  result  That  is  why  the  sentence  img  =  -img  is  used  in  the 
program  calling  the  routines. 

In  order  for  the  set  of  samples  of  a  signal  to  correctly  represent  that  signal,  the 
sampling  theorem  has  to  be  satisfied.  Since  an  upper  limit  of  integration  or  summation  xmax 
has  been  assumed,  the  sampling  frequency  can  not  be  less  than  2xmax.  Suppose  we  want 
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L  points  per  cycle  for  the  incident  burst,  then  the  sampling  theorem  states  the  following 
condition 

Lx02  2xmax.  (A6.1) 

Throughout  the  report,  xmax  ■  400  and  Ax  =  0.05  have  been  chosen  except  for  the 
aluminum  and  aluminum/fluid  shell  cases  where  xmax  =«  480  and  Ax  *  0.1.  When  xq  <  40, 
we  can  either  let  L  be  greater  or  reduce  the  value  of  xmax  a  little  bit  We  decided  to  adopt  the 
former  solution,  namely,  changing  L  while  keeping  xmax  =  400.  Examples  are  Fig.6d 
where  L  »  40  and  Fig  5c(a),  1 1  a,  1 1  b,B  1  (a)  and  B 1  (b)  where  L  *  25. 

As  stated  in  Chapter  4,  large  value  of  N  could  be  a  calculation  problem.  A  10- 
cylinder  account  on  the  CMS  main  frame  can  only  support  the  calculation  of  N-34000. 
When  computing  the  backscattering  from  aluminum  or  aluminum/fluid  shell,  the  carrier 
frequencies  we  are  interested  in  are  Xq  *  302  and  340.  Since  we  have  N  =  Lx^/Ax,  L  has 
to  be  decreased  and  Ax  has  to  be  increased  in  order  to  keep  N  under  this  limit  Therefore 
Ax  »  0.1  and  L  =  10  have  been  chosen,  which  also  assures  the  sampling  theorem 
condition  for  xmax  *  480. 

It  should  be  pointed  out  that  the  upper  limit  of  the  summation  in  Eq.(4.15)  should 
be  the  same  as  xmax.  In  other  words  AxN/2  »  xmax,  or  Lx0  *  2xmax.  Fortunately,  the 
components  above  xmax  for  all  the  bursts  considered  in  this  report  is  small  enough  that  as 
far  as  the  physics  in  the  scattering  process  is  concerned,  the  inequality  (A6.1)  can  be 
regarded  as  the  proper  condition  which  has  been  assured  throughout  the  work. 


imperfection  is  mainly  at  the  two  edges. 


B.  Additional  Data  and  Graphs 


Bl.  Mid-frequency  enhancement  of  the  backscattering  for  4-cycle  bunt 
B2.  Resonance  backscattering  from  aluminum/fluid  shell:  1 -cycle 
B3.  Resonance  backscattering  from  aluminum/fluid  shell:  4-cycie 
B4.  Resonance  backscattering  from  SS440c/fluid  shell:  1 -cycle 
B 5.  Resonance  backscattering  from  SS440c/fluid  shell:  4-cycle 
B6.  Resonance  backscattering  from  real  aluminum  shell:  1 -cycle 
B7.  Resonance  backscattering  from  real  aluminum  shell:  4-cycle 
B8.  Resonance  backscattering  from  real  SS440c  shell:  1 -cycle 
B9.  Resonance  backscattering  from  real  SS440c  shell:  4-cycle 


BACKSCATTERING  FROM  SPHERICAL  SHELL 


BACKSCATTERING  FROM  SPHERICAL  SHELL 


BACKSCATTERING  FROM  SPHERI 


Fig.  Bl(c).  Computational  backscattcring  from  a  2.5%-thick  SS304  shell  for  a  4- 
cycle  sine  burst  with  xq  =  42.0. 


Fig.  Bl'd).  Computational  backscattering  from  a  2.5%-thick  SS304  shell  for  a  4- 
cycle  sine  burst  with  xo  =  46.0. 


Fig.  Bl(e).  Computational  backscattering  from  a  2.5%-thick  SS304  shell  for  a  4- 
cyclc  sine  burst  with  xq  =  50.0. 


Fig.  Bl(f).  Computational  backscattering  from  a  2.5%-thick  SS304  shell  for  a  4- 
cycle  sine  burst  with  Xq  =  54.0. 


Fig.  Bl(g).  Computational  backscattering  from  a  2.5%-thick  SS304  shell  for  a  4- 
cycle  sine  burst  with  x0  =  58.0. 


Fig.B2(b).  Computational  backscattering  from  a  4%-thick  aluminum/fluid  shell  for 
a  1 -cycle  sine  burst  with  xo  =  340. 


Fig.B3(a).  Computational  backscattering  from  a  4%-thick  aluminum/fluid  shell  for 
a  4-cycle  sine  burst  with  xq  =  302. 


Fig.B3(b).  Computational  backscattering  from  a  4%-thick  aluminum/fluid  shell  for 
a  4-cycle  sire  burst  with  xq  =  340. 


BACKSCATTERING  FROM  SPHERICAL  SHELL 
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Fig.B4(a).  Computational  backscattering  from  a  16.2%-thick  SS440C/fluid  shell 
for  a  1 -cycle  sine  burst  with  x0  =  71.0. 


BACKSCATTERING  FROM  SPHERICAL  SHELL 
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Fig.B4(b).  Computational  backscattering  from  a  16.2%-thick  SS440C/f1uid  shell 
for  a  1-cycle  sine  burst  with  xo  =  76.8. 


BACKSCATTERING  FROM  SPHERICAL  SHELL 


Fig.B5(a).  Computational  backscattering  from  a  16.2%-thick  SS440C/fluid  shell 
for  a  4-cycle  sine  burst  with  xo  =  71.0. 


BACKSCA7TERING  FROM  SPHERICAL  SHELL 
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Fig.DS(b).  Computational  backscattering  from  a  16.2%-thick  SS440C/fluid  shell 
for  a  4-cycle  sine  burst  with  xq  =  76.8. 


Fig.B6(a).  Computational  backscattering  from  a  4%-thick  real  aluminum  shell  for 
a  1 -cycle  sine  burst  with  xq  =  302. 


Fig.B6(b).  Computational  backscattcring  from  a  4%-thick  real  aluminum  shell  for 
a  1 -cycle  sine  burst  with  xq  =  340. 


ig.B7(a).  Computational  baclcscattering  from  a  4%-thick  real  aluminum  shell  for 
4-cycle  sine  burst  with  xr.  -  302. 


Fig.B7(b).  Computational  backscattering  from  a  4%-thick  real  aluminum  shell  for 
a  4-cycle  sine  burst  with  xQ  =  340. 


Fig.B8(b).  Computational  backscattering  from  a  16.2%-thick  real  SS440C  shell 
for  a  1-cycle  sine  burst  with  xo  =  76.8. 


Fig.B9(a).  Computational  backscattering  from  a  16.2%-thick  real  SS440C  shell 
for  a  4-cyclc  sine  burst  with  x©  =  71.0. 


BACKSCATTERING  FROM  SPHERICAL  SHELL 


Fig.B9(b).  Computational  backscattering  from  a  16.2%-thick  teal  SS440C  shell 
for  a  4-cycle  sine  burst  with  Xq  “  76.8. 
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C.  Impulse  Response  and  Form  Function 

Eq.(3.7)  states  that  p’(t)  and  F(x)  are  a  Fourier  Transform  pair  when  G(x)  a  1,  or 
in  other  words,  when  incident  wave  is  a  Dirac  delta  function  5(x).  Therefore  we  would 
expect  to  obtain  the  impulse  response  of  the  scatterer  by  taking  the  inverse  Fourier 
Transform  of  the  form  function.  Considering  the  time  reference  point,  we  let  the  incident 
wave  stardng  at  x  =  3  be  a  half  cycle  sine  wave  with  xq  =  200  and  L  =  4.  It  is  equivalent  to 
a  delta  function  with  nonzero  value  1  only  at  t  =  3.  Its  spectrum  was  checked  to  be  unity 
for  all  x  value.  When  taking  the  inverse  Fourier  Transform,  the  cut-off  frequency  is  chosen 
to  be  400.  Computed  backscattering,  or  the  impulse  response  of  a  2.5%-thick  SS304  shell 
is  shown  in  Fig.Cl.  We  clearly  see  a  low  frequency  response  starting  at  t ■  5  which  is  the 
manifestation  of  existence  of  a  mid-frequency  enhancement.  The  measured  interval 
between  the  center  of  the  response,  assumed  to  be  right  at  the  peak  indicated  by  an  arrow, 
and  the  specular  reflection  is  4.44.  The  excellent  agreement  with  the  ray  model  calculation 
value  (4.467)  strongly  suggests  the  validity  of  the  coincident  ray  model. 

Fig.  5c(a),  5c(b)  and  5c(c)  suggest  a  very  low  frequency  enhancement.  The 
impulse  response  of  a  16.2%-thick  SS440C  shell,  shown  as  Fig.C2,  was  computed  by 
choosing  a  half-cycle  sine  burst  with  xo  =  200  and  L  =  4.  A  low  frequency  response  is 
evident,  which  shows  that  there  does  exist  a  low  frequency  enhancement  for  such  a  thick 
shell. 
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D.  Additional  Evidence  of  Mid-Frequency  Enhancement: 

Arrival  Time  Check  With  Ray  Model  Calculation 

The  coincident  ray  model  for  the  mid-frequency  enhancement  has  predicted  that 
relative  to  the  specular  reflection  the  arrival  time  of  the  echo  due  to  the  mid-frequency 
enhancement  can  be  described  by  td  *  2  +  Jtc/cg/,  The  calculated  time  interval  between  the 
center  of  the  specular  reflection  and  the  center  of  the  coincident  echo  by  ray  model  is  4.467 
for  2.5%-thick  SS304  shell  and  xo  -  46.0.  Fig. lid  shows  that  at  x<j  *  46.0,  the  echo  due 
to  die  mid-frequency  enhancement  is  quite  symmetric  therefore  it  is  easy  to  locate  the  center 
of  the  echo.  Actually  it  can  be  determined  as  t  ■  6.768  by  looking  at  the  Table  4.4  and  the 
center  of  the  specular  reflection  can  be  calculated  as  t  =  1  +  10  x  2n/x0  *  2.366.  The  time 
interval  between  them  is  therefore  4.402  which  is  pretty  much  the  same  as  the  ray  model 
prediction.  This  simple  check  gives  additional  evidence  of  the  validity  of  the  coincident  ray 
model. 


t 


Fig.Dl.  The  arrival  time  of 
various  surface  waves  for 
2.5%-thick  SS304  shell.  The 
solid  curve  is  for  /  =  ao,  the 
short-dashed  curve  is  for  /  = 
so,  and  the  short-long-dashed 
curve  is  for  /  =  ao'. 
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E.  More  on  Prompt  Radiation  Effect 

We  have  seen  in  section  5.5  that  the  prompt  radiation  can  be  out  of  phase  with  the 
specular  reflection,  which  makes  it  difficult  to  infer  the  amplitude  of  the  prompt  radiation 
wave  for  the  case  of  the  backscattering  from  a  SS440C  shell  (Fig.l8a).  In  section  5.4  we 
have  assumed  a  fluid  shell  to  remove  the  effects  due  to  the  elastic  surface  waves  and 
transverse  waves.  Here  we  want  to  single  out  the  prompt  radiation  from  the  specular 
reflection  background,  which  can  be  achieved  by  subtracting  the  backscattering  by  a  fluid 
shell  from  that  by  an  elastic  shell.  The  results  for  aluminum  and  SS440C  shell  cases  have 
been  shown  as  Fig.E.  The  prompt  radiation  can  be  clearly  seen  for  both  aluminum  (xq  = 
302)  and  SS440C  (xo  *  71.0)  cases.  Comparison  of  Figs.  El  and  E2  for  the  aluminum 
case  shows  clear  evidence  of  a  substantial  prompt  contribution  for  xo  =*  302  and  a  much 
weaker  prompt  contribution  for  x0  =  340.  It  is  interesting  to  notice  that  Figs.  El  and  E2 
display  a  qualitative  similarity  with  the  burst  response  of  a  damped  harmonic  oscillator  near 
and  away  from  resonance  respectively.  For  a  review  of  the  burst  response  of  such  systems 
see  e.g.  Pippard23. 

The  arrival  of  other  surface  wave  contributions  makes  the  steady-state  prompt 
radiation  contribution  hard  to  be  determined  for  the  SS440C  case.  Shown  in  Fig.E5  is  the 
case  of  80-cycle  incident  burst.  If  the  backscattering  reaches  its  steady  state  right  before  the 
burst  nuns  off,  which  is  a  suggestion  of  the  figure  although  it  can  not  be  fully  trusted,  then 
the  measured  amplitude  value  is  2.0.  Neglecting  the  relatively  weak  contribution  of  other 
surface  waves  in  this  time  interval,  the  value  2.0  becomes  the  estimate  of  the  prompt 
radiation  contribution.  This  value  is  remarkably  close  to  the  value  of  1.95  at  Xo  =  71.0 
implicit  in  Fig.7  of  Ref.20.  That  is  the  magnitude  of  the  error  in  the  steady-state  ray 
synthesis,  which  neglects  prompt  radiation  contribution.  Thus  1.95  and  2.0  are  consistent 
estimates  of  the  prompt  contribution  for  the  SS440C  case. 
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Fig.E5.  The  difference  between  the  backscattering  by  a  16.2%-thick  SS440C  shell 
and  that  by  a  16.2%-thick  SS440C/fluid  shell  for  a  80-cycle  sine  burst  with  Xn  = 
71.0. 
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